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A mixed Monte Carlo and PDE variance reduction method 
for foreign exchange options under the Heston-CIR model 

ANDREI COZMA* CHRISTOPH REISINGER* 


Abstract 

In this paper, the valuation of European and path-dependent options in foreign exchange 
(EX) markets is considered when the currency exchange rate evolves according to the 
Heston model combined with the Cox-Ingersoll-Ross dynamics for the stochastic domes¬ 
tic and foreign short interest rates. The mixed Monte Carlo/PDE method requires that 
we simulate only the paths of the squared volatility and the two interest rates, while 
an “inner” Black-Scholes-type expectation is evaluated by means of a PDE. This can 
lead to a substantial variance reduction and complexity improvements under certain 
circumstances depending on the contract and the model parameters. In this work, we 
establish the uniform boundedness of moments of the exchange rate process and its ap¬ 
proximation, and prove strong convergence in (p > 1) of the latter. Then, we carry 
out a variance reduction analysis and obtain accurate approximations for quantities of 
interest. All theoretical contributions can be extended to multi-factor short rates in 
a straightforward manner. Einally, we illustrate the efficiency of the method for the 
four-factor Heston-CIR model through a detailed quantitative assessment. 

Keywords: conditional Monte Carlo, mixed Monte Carlo/PDE, stochastic volatility, 
stochastic interest rates, variance reduction, strong convergence. 


1 Introduction 


In FX markets, option pricing with stochastic volatility and stochastic interest rates has seen 
a large amount of interest in the last few years ( Grzelak and Oosterlee|20fT ; Van Haastrecht 


and Pelsserl|2011 ; Ahlip and Rntkowski||2013 ), leading to the extension of the Heston (1993) 


two-factor stochastic volatility model and the Schobel-Zhu ( Schobel and Zhu|1999 ) model to 
currency derivatives. Although appealing due to its simplicity, assuming constant interest 
rates is inappropriate for long-dated FX products, and the effect of interest rate volatility 
can even outweigh that of FX rate volatility for long maturities, a fact confirmed by empirical 
results (Van Haastrecht et al. 2009). Here, the spot FX rate is defined as the number of 
units of domestic currency per one unit of foreign currency. 
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In this paper, we consider the four-factor Heston-CIR model proposed and examined 


in Ahlip and Rutkowski (2013) when the volatility and the exchange rate dynamics are 
correlated, whereas the domestic and foreign interest rates are pairwise independent and 
also independent of the exchange rate and the volatility. Our motivation comes from the 
fact that the square root (CIR) process ( Cox et aT]|1985 ) for the variance and interest rates 
is widely used in the industry due to its desirable properties, such as mean-reversion and 
non-negativity. Under these restrictive assumptions on the independence of the Brownian 
drivers, the authors argue that the model is affine and derive a semi-analytical formula for 
the European call option price. The importance of a non-zero correlation between the FX 
rate and the interest rate(s) is recognized in Hunter (2005). However, any other non-zero 
correlations give rise to a non-affine model, in which case we lose analytical tractability. 

Therefore, we see ourselves forced to turn to numerical algorithms, and the mixed Monte 
Carlo/PDE solver ( Loeper and Pironneau||2009 ) is a good alternative to the classical Monte 
Carlo and finite difference methods. Monte Carlo simulation methods (Glasserman 2003) 
can handle path-dependent features easily and scale linearly with the dimension, however 
a considerable number of simulations is typically required for a good accuracy. Conversely, 
finite difference methods incorporate early exercise features easily and provide a fast con¬ 
vergence for low-dimensional problems (up to three dimensions), but become intractable 
as the dimensionality increases. This makes Monte Carlo methods more attractive for the 
four-factor Heston-CIR model, and the relatively slow convergence in the number 

of simulations M raises the question of finding an efficient variance reduction technique. 

The idea behind the mixed Monte Carlo/PDE method is to write the option values as 
nested conditional expectations. Then, the innermost expectation is evaluated analytically, 
in case of European-style options, or using finite differences, and the outer expectation by 
simulation. The earliest related published work belongs to Hull and White (1987). The au¬ 
thors consider a European call option under a two-dimensional stochastic volatility model 
with uncorrelated asset price and volatility and prove that, conditional on the integral of 
the variance process, the asset price is lognormally distributed and hence the option price 


can be expressed as a Black-Scholes price. Willard (1997) extends the analysis of Hull and 


White (1987) to path-independent options under stochastic volatility and instantaneously 
correlated factors, and uses the smoothness of the “conditional price” to calculate price 
sensitivities (the Greeks). The author also employs quasi-Monte Carlo (low-discrepancy) 
methods to further reduce the variance of price estimates, but with no effect on the dis¬ 
cretization bias. The mixed Monte Carlo/PDE method develops this conditioning technique 
- known as conditional Monte Carlo - and allows the combined use of Monte Carlo and 
finite difference methods for the valuation of path-dependent contracts. 

The utility of the method can be easily recognised when pricing European-style options 
under the Heston-CIR model. Conditioning on the entire paths of the squared volatility and 
the domestic and foreign interest rates, the dynamics of the exchange rate are governed by a 
geometric Brownian motion with time-dependent drift and diffusion coefficients. Combined 
with the existence of a closed-form solution for the conditional option price, the algorithm 
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reduces the variance in Monte Carlo simulations - by eliminating a source of noise - and the 
dimension of the problem, from four to three. The mixed Monte Carlo/PDE method has 
been the object of several numerical studies in the past few years, e.g., |Lipp et si] ( |2013[ ), 


Ang (2013), McGhee (2014), Dang et al. (2015), and various extensions to the original idea 


have been considered. For instance, Dang et al. (2015) use the Hull-White type dynamics 
of the interest rates and condition on the variance path to find a closed-form solution for 
the conditional price of a European option. A few of these references examine convergence 
properties of the algorithm by heuristic arguments, and none of them take into account the 
error arising from the discretization of the variance and interest rate processes. 

To the best of our knowledge, the convergence of the mixed Monte Carlo/PDE method 
has not yet been established, and even the literature on Monte Carlo methods under stochas¬ 
tic volatility is scarce. Higham and Mao (2005) considered an Euler simulation of the Heston 
model with a reflection fix and proved strong convergence of the stopped approximation pro¬ 
cess, as well as for a European put and an up-and-out call, by using the boundedness of 
payoffs. Cozma and Reisinger (2015a) extended these results to derivatives with unbounded 
payoffs, stochastic-local volatility, stochastic interest rates, and exotic payoffs. Their dis¬ 
cretization scheme coincides with the one considered in this paper at the discrete time 
points, but necessarily differs in the continuous-time interpolation. Some results can be 
utilized from the earlier work, although stronger conditions on the model parameters were 
needed there and therefore the results here can also be seen as an improvement for the stan¬ 
dard (i.e., non-mixed) scheme. In particular, the new results are always guaranteed for zero 
or negative correlation in the Heston model. The main conceptual difference, however, is 
the new continuous-time interpolation using conditional drifts, which is crucial for deriving 
the conditional PDFs and leads to different technical challenges. Also, for path-dependent 
options, neither the original scheme nor the analysis applies to the present work. 

In this paper, we study convergence properties of the mixed Monte Carlo/PDE method 
with the full truncation Euler (FTE) discretization for the squared volatility and the two 
interest rates, and demonstrate the efficiency of the method for a European call and an 
up-and-out put option. We prefer the full truncation scheme (Lord et al. 2010) because it 
preserves positivity, is easy to implement and is found empirically to produce the smallest 
bias among all Euler schemes. An interesting alternative to the FTE scheme would be the 
backward Euler-Maruyama (BEM) scheme (Neuenkirch and Szpruch 2014). However, the 
quanto correction term in the dynamics of the foreign interest rate would lead to technical 
challenges in the convergence analysis. The major contributions of this paper are as follows. 


• We establish the uniform boundedness of moments of the four-dimensional process and 
its approximation, and prove strong convergence in (p ^ 1) of fho discretization 
scheme. Then, we deduce the convergence of mixed Monte Carlo/PDE estimators for 
computing option prices and discuss possible extensions to higher-dimensional models. 

• We carry out a thorough theoretical variance reduction analysis of the mixed Monte 
Carlo/PDE method and employ standard Monte Carlo - with the log-Euler discretiza- 
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tion - as the reference method, noting that the analysis applies to general interest rate 
dynamics. In particular, we investigate how different values of the underlying model 
parameters affect the variance of mixed estimators. 

• We perform a series of numerical experiments and demonstrate the convergence of the 
mixed Monte Carlo/PDE method under the four-factor FX model for two financial 
derivatives: a European call and an up-and-out put option. In addition, we establish 
the efficiency of the method by comparison with alternative numerical schemes and 
examine the sensitivity of the variance reduction factor to changes in the parameters. 

The remainder of this paper is structured as follows. In Sectionwe introduce the four- 
factor FX model, define the mixed simulation scheme and describe the pricing algorithm. 
In Section we prove strong convergence of the exchange rate approximations and then 
discuss some extensions. Detailed proofs of some technical results are given in the Appendix. 
In Section]^ we carry out a variance reduction analysis of the mixed Monte Carlo/PDE 
method for a European option. Various numerical experiments are presented and discussed 
in Section Finally, Section summarizes the results and outlines possible future work. 


2 Preliminaries 

2.1 The four-factor model 

We have in mind a model in an FX market, for the spot FX rate 5, the variance of the FX 
rate the domestic short interest rate and the foreign short interest rate . Unless oth¬ 
erwise stated, in this paper, the subscripts and superscripts “d” and are used to indicate 
domestic and foreign, respectively. Consider a filtered probability space (D, Q) 

and suppose that the dynamics of the underlying processes are governed by the system of 
stochastic differential equations (SDEs) below under the domestic risk-neutral measure Q: 

'dSt = {rf - r{)Stdt + ^/^HStdW^ 
dvt — k(^0 — vt)dt + ^^Jvt dWi 

< ^ ( 2 . 1 ) 

dvi = kdiOd - rf)dt + Uyrt dWf 

Jr{ = {kfOf - kfr{ - psf^f\fv^t)dt + ^f^dw /, 

where VF^, } are correlated standard Brownian motions (BMs) under the risk- 

neutral measure with constant correlation matrix X. We consider a full correlation structure 
between the four Brownian drivers, which reflects movements in the financial markets more 
accurately and allows for better calibrations. Then, we decouple {VF^, VF-^, VF^, VF^} and 
express them as linear combinations of independent Brownian motions {VF^, VF^, VF^, W^}. 
Hence, define the two vectors W = [VF", VF^, VF^, and W = [VF\ VF^, VF^ As 
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S is symmetric positive definite, a standard Cholesky factorisation gives rise to an upper 
triangular matrix of coefficients A = ((^ij)i<ci j< 4 ^ satisfying S = ^4^4^, which is given below. 
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0 

022 

023 

024 

Psd 

Pdf 

1 

pvd 


0 

0 

O33 

O34 

Psv 

Pvf 

Pvd 
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0 

1 


This decomposition implies that we can choose W so that W = AW. We can determine the 
matrix of coefficients by solving a system of ten equations. Assuming ^ ±1, we find: 


<^14 — Psv^ ^24 — Pvf^ ^34 — Pvd^ ^33 — (l Pvd) ^ ’ ^13 — (Psd PsvPvd) (l Pvd) ^ ’ 

<^23 — {pdf ~ PvfPvd) (l ~ Pvd) ^ ’ ^22 ~ (l ~ Pdf ~ Pvf ~ Pvd ‘^Pdf Pvf Pvd) ^ (l ~ pyd) ^ , 


whereas an and ai 2 can be found in a similar fashion. 

The quanto correction term in the drift of the foreign interest rate in comes from 
changing from the foreign to the domestic risk-neutral measure ( Clark|2011[ ). Alternatively, 
we can also think of (2.1) as a model in an equity market with asset price process interest 
rate and dividend yield r-^, in which case the quanto drift adjustment term vanishes. 


2.2 The mixed simulation scheme 


First, 

Euler 


we discretize the variance and the two interest rate processes using the full truncation 
(FTE) scheme of Lord et al. (2010). Consider the square root process 


dyt = kyiOy - yt)dt + dW^. 


(2.3) 


For a time interval [0,T], consider a uniform grid: 5t — T/N^ tn = Vn G {0,1,..., A^}. 
We introduce the time discrete auxiliary process 


yt„+i = ytn + ky {By - yf)5t + 4 , (2.4) 

where = max(0,p) and — W^^^^ — Wf^^ and the time continuous interpolation 


Vt — ytn + ^yi^y ytn){^ ^y\/ytn{^t ^tn) ^ ^ [tuAn+l)^ (2-5) 


as suggested in Higham and Mao (2005). Moreover, we define the non-negative processes 


and 


Yt = yt 

(2.6) 

Yt = yt, 

(2.7) 
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whenever t G Let V and be the FTE discretizations - as defined in (2.7) - of 

the variance and the domestic interest rate, respectively. Taking into account the presence 
of the quanto correction term in the drift of the foreign interest rate, we similarly define 




*/«/- (t - tn) + - w'/j, ( 2 . 8 ) 


as well as 

H = {H)^ (2-9) 

and 

H ( 2 - 10 ) 

whenever t G Next, we define 5, the continuous-time approximation of 5, as the 

solution to the following SDE: 


dSt = iJ.tStdt + aiiyVtStdWt, 


Ait = rf - rf - t (1 - a^i) Vt + 


i=2 


7= 

6t 


( 2 . 11 ) 


where 6W^ — — kE/^, Vt G [t^, t^+i), hence /it is piecewise constant. For convenience, 

we introduce the actual and the approximated log-processes, x = log S and X = log S. 

Conditioning on the trajectories of {kE-^, j = 2,3,4}, i.e., on the complete knowledge 
of the paths of the variance and interest rates, S evolves like a geometric Brownian motion 
with time-dependent drift and diffusion coefficients. It follows from Ito’s formula that 

5T = 5oexp|y {ri-fi--Vu)du + Y,aij j ^JVu^du + an j X-^dWl 

j_2 


However, we know the conditional probability law of the stochastic integral on the right- 
hand side to be that of a normal random variable, namely 



law 



Vudu • Z^ 


where Z ^ A/” (0,1), so we can think of St as a function of Z. Therefore, we can express it 
as 


St = 5oexp | (^r - g - + cr\/rz| , 


( 2 . 12 ) 


where 


r = 


1 

T 



rfi^du 


1 

N 


N-l 




i=0 


a 



N-l 


^11 

N 




i=0 
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J =2 
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Conditional on the trajectories of {W^^j = 2,3,4}, ( |2.12 ) has the same law as a terminal 
asset price that evolves as a geometric Brownian motion with interest rate r, continuous 
dividend yield q and volatility a, all constant. Then the arbitrage-free price at time t = 0 
of a European-style option with payoff /{St) is the discounted expectation under the risk- 
neutral measure, which can be approximated using the mixed Monte Carlo/PDE method 
by 


U = E 


e-fo 


= E 


E 


o-rT 


f(ST)\a 


fAv 

T 


(2.13) 


where 0 < t < T} is the natural filtration generated by the independent Brownian 

motions VE^}, i.e., generated by the processes i;, r^, as observed until time T. 

The second equality in ( |2.13 ) comes from the “tower property” of conditional expectations. 
Unless otherwise stated, all expectations are under Q. Let the approximate conditional 


option price be the inner expectation in (2.13), which is analytically tractable for European 
contracts. 


E 


1 ^oo 

,-rrj(^^)| ^g-rT / f(^ST{z))(f>{z)dz, 

J —oo 


where (j) and $ are the standard normal PDF and CDF, respectively. The conditional prices 
of some popular financial instruments are given below. 


{ European options: ^{/jdi) — ^{i/>d 2 ) 

Cash-or-nothing options: ^{/jd 2 ) 

Asset-or-nothing options: 

where = 1 for a call and = — 1 for a put, whereas 

\og{S^/K) + {r-q±cj^/2)T 


(2.14) 


G?1,2 = 


a 


Vf 


(2.15) 


The approximate option price, i.e., the outer expectation in (2.13), is estimated by a Monte 


Carlo average over a sufficiently large number of discrete trajectories of {VE^, VE^, VE^}. 

There are many other derivatives that admit a closed-form solution for the conditional 
price, like the power option, the chooser option or the forward-start option. However, for 
most path-dependent derivatives, we need to use a different approach in order to compute 
the conditional price, in which case we will rely on finite difference methods (see Section]^. 
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We conclude this section with a discussion on our choice of conditioning. For European 
option pricing, an analytical formula for the inner expectation is available only when condi¬ 
tioning on all three factors, which results in a dimension reduction of the problem by one. 
The high-dimensionality of the Heston-CIR model makes this the natural choice. For path- 
dependent option pricing, due to the quanto correction term in the drift of the foreign rate, 
V and are coupled and hence we cannot condition on alone. Moreover, the variance of 
Monte Carlo estimators due to the short rates is typically much lower than the variance due 
to the instantaneous squared volatility. Hence, we could alternatively condition on and 
solve a three-dimensional PDE for the inner expectation. On one hand, we could reach the 
same level of accuracy with fewer Monte Carlo sample paths than when simulating v and 
as well. On the other hand, the computational effort grows linearly with the dimension 
for Monte Carlo methods and exponentially for finite difference methods. Hence, we believe 
that conditioning on all three factors is more efficient. However, if the exchange rate and 
the foreign interest rate dynamics are independent, the quanto correction term vanishes. In 
this case, conditioning on the two short rates would be an interesting alternative. 


3 Convergence analysis 

Even though weak convergence is very important in financial mathematics when estimating 
expectations of payoffs, strong convergence may be required for complex path-dependent 
derivatives and plays a key role in multilevel Monte Carlo methods ( Giles|2008 ). In this sec¬ 


tion, we prove the strong convergence of the approximation scheme defined in (2.11). Hence, 


we first examine exponential integrability properties of the square root process and its dis¬ 
cretization, and then the finiteness of moments of order higher than one of the exchange 
rate process and its approximation. 


Let y be the square root process defined in (2.3) and let Y be the piecewise constant FTE 


interpolant from (2.7). The exponential integrability of functionals of the two processes was 


already discussed in Propositions 3.2 and 3.6 in Cozma and Reisinger (2015b). However, we 


need to adjust their second result for our approximation scheme (2.11) in order to establish 
the convergence. 


Lemma 3.1. Let A,/x G M he given, A = A -h , and define the stoehastie proeess 


Qt = exp I ^ Yudu + y j 


6t 


duY Vt G [0,T]. 


(3.1) 


//A < 0 and T >0, or otherwise, i/A > 0 and T then there exists 77 > 0 sueh that 


sup sup E [Qt] < ex), 
ste{0,r]) telo,T] 


(3.2) 


where T* is given below: 











1. Ifky<iy{^Ji+V^:^), 


2 . Ifky > ^y(M + Vo. 5 A), 


Proof. See Appendix [A| 


_ 


^^(/i + \/M) - ky 


rp^ _ 2(fcy /i(^y) 


(3.3) 


(3.4) 


3.1 Moment bounds 

For many stochastic volatility models, moments of order higher than one can explode in 
finite time ( Andersen and Piterbarg|[2007 ). This can cause significant problems in practice, 
for instance when computing the arbitrage-free price of an option whose payoff function 
has super-linear growth. The same troublesome behaviour can be observed for the Euler- 
Maruyama approximation of some SDEs with super-linearly growing drift or diffusion co¬ 
efficients, where moments diverge in finite time ( Hntzenthaler and Jentz^|2015 ). Next, we 
prove the boundedness of moments of the exchange rate process and its approximation. 

At this point, we assume that ^ ±1 and that ai 3 is non-zero, i.e., psd ^ PsvPvd- 

Proposition 3.2. For a > 1, define the two quantities 

1 


qo{a) = 




13 


qi{a) = qo{a) lp,„<o + min qo{a), 


2 apsv^k + + 0^2) — 

- [ 2 apsv^k + + 0^2) - 

k 


(^Psv^ 


1 


Psv ^ 0 ' 


If the following eonditions on the model parameters are satisfied, 

1 ^ r-p -7 k‘^ aqiia) 

k > apsvf, + V <^(<^ ~ 1 ) C > 


qi{a)-l' 

then there exists ai > a such that for all u € [ 1 , ai), 

sup E [S”^] < 00 . 
te[o,T] 

Proof. See Appendix |B| 

Next, assume that 013 and 014 are not simultaneously zero, i.e., 7 ^ 0. 


(3.5) 

(3.6) 


(3.7) 


(3.8) 
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Proposition 3.3. For a > 1, define 


q 2 {a) = + al^)/A - Tae/^‘' + Ta^^{al^ + a\^)k 

- [ap,^^ + Ta^i^{a\^ + al^)/A - Taf/A\ \ 2 , 2 n • ( 3 - 9 ) 

) la -h 

If the following eonditions on the model parameters are satisfied, 

k> aps^^+^a{a-l)T^‘^ and ^ ^ «g2(«) ^ 

4 rq g2(a) - 1 

then there exists a 2 > a sueh that for all cj G [1, 0 ^ 2 ); "we ean find 77 ^ > 0 50 that 

sup sup E [S^] < 00 . (3.11) 

Ste{0,r]u^) tG[0,T] 

Proof. See Appendix \C\ ■ 

Since the most popular FX and equity contracts grow at most linearly in FX and asset 
prices, and their valuation requires the computation of the expected discounted payoff under 
the risk-neutral measure, it is useful to study finiteness of moments under discounting. Let 
R be the discounted exchange rate process, 

= So exp I - y (rl + ^Vufiu +J ^/vZdwfj, (3.12) 

and let R be its continuous-time approximation, 

Rt = So exp I - + i Vu) du + an dW^ + ^ ay y^ ^ • (3-13) 

J —2 

Proposition 3.4. Let a>l. IfT < T*, there exists ai > a sueh that for all u G [l,ai), 


sup E [i?^] < oo. 
te[o,T] 


If a > 1 and T >T*, then 

E [Rif] = oo. 

If a = 1, then T* = oo, whereas if a > 1, then T* is given below: 
1. If k < apsvi - \/a{a - 1)^, 

T* = G_ log ( 

v{a) ^\aps^^ - k - v{a)) ' 


(3.14) 


(3.15) 


(3.16) 
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where 


iy{a) = V {apsv^ - k)^ - a{a - 1)^2 _ 
2. If k = apsvi - \/a(a - 1)^, 


_ 


^Psv^ ^ 

3. If apsv^ - y/a{a - 1) ^ <k < apsv^ + \/ a{a - 1) 


(3.17) 


_ 


v{a) 


f otpgvf, k 


TV 

-arctan i — ,, , 

2 V ^>(«) 


(3.18) 


where 


P{a) = a/ a{a - 1)^2 _ (ap^^^ _ fc)2. 

4- If k> apsvi + ^/a{a - 1 )^, 


T* = oo. 


(3.19) 


Proof. See Appendix [D| ■ 

When the domestic and the foreign interest rates are constant, Proposition |3.4| examines 
moment boundedness in the Heston model. It is an extension of Proposition 3.1 in | Andersen 


and PiterbargI ( |2007 ) from bounds on moments of order a > 1 to bounds on all moments of 


order u G [a, ai), for some ai > a > 1. We use this result to prove the strong convergence 
of the discretized discounted spot FX rate process. 


Proposition 3.5. Let a > 1. If T < T*, there exists a 2 > a and > 0 sueh that for all 
id G [1, 0 ^ 2 ); 

sup sup E \B!f] < 00 , (3.20) 

ste{o,r]u;) telo.T] 

where T* is given below: 

1 . Ifk< apsv^ + \ 

r* =- ^ -. (3.21) 

apsvi + V 


2. Ifk> apsvi + \ 




OO 

4:{k - apsvQ 
a{a - 1)^2 


if a = 1 
if a > 1. 


Proof. See Appendix H 


(3.22) 


11 






















To the best of our knowledge, the boundedness of moments of discretization schemes for 


the Heston model and extensions thereof had not been established until recently (Cozma 


[and Reisinger|2015a ) - a fact that is also mentioned in Kloeden and Neuenkirch (2012) - and 
Proposition |3.5| is only the second to address this issue. For the Heston model, Proposition 


|3.5| can be seen as an improvement of Proposition 3.9 in Cozma and Reisinger (2015a) due 
to the sharper conditions on the critical time. 

3.2 The four-dimensional system 

The strong mean square convergence of the discretized variance and domestic interest rate 


processes was established in Proposition 3.5 in Cozma and Reisinger (2015a). First of all. 


we prove an equivalent result for the foreign interest rate. 

Proposition 3.6. If2kf9f > then the proeess eonverges strongly in L^, i.e., 


-f 

lim sup E 
st^o tG[o,r] 


I / -/|2' 

In - n \ 


= 0 . 


(3.23) 


Proof. See Appendix]^ 


As will become clear from the proof, the Feller condition 2kf0f > which ensures that 
the process r^ does not hit zero, allows us to control the potential growth of the absolute 
difference between the original and the discretized processes that comes from the sublinear 
correction term in the drift. 


Second, we consider the logarithm of the process from (2.11) and examine its convergence 


properties. The formulae of the log-process, x, and its approximation, X, are given below. 


Xt = Xo + 


(rt-rl-^Vu^du + ^dW^, 


rt 2 _ \ f 

Xt = XQ+ / [f^-fl--Vujdu + an / yK dW^ + ^ aij / 

^ 0 ^ 0 j 2 ^ 0 


(3.24) 


Vu du. (3.25) 

/o ' 


Proposition 3.7. If2kf6f > then the log-proeess eonverges uniformly in , i.e., 


lim E 


sup \xt — Xt 


= 0 . 


(3.26) 


lte[o,T] 

Proof. See Appendix ■ 

Third, we prove convergence of the discretized spot FX rate process. 

Proposition 3.8. If2kf9f > then the proeess S eonverges uniformly in probability, i.e.. 


lim P ( sup > e ) = 0, Ve > 0. 

St —>-0 \tG[0,T'] 


(3.27) 
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Proof. See Appendix [H| 

Theorem 3.9. Let a > 1 and assume that the following eonditions are satisfied: 
k > apsv^ + max | ~ I, 


k‘^ aqi{a) 2kd aq 2 {a) 


> 


TT > —-- and 2kf9f > 

2 g2(a)-l ^ 


2 {J 9l(tt)-l’ Til 

Then the proeess eonverges strongly in in the sense that 


lim sup E 

St —>-0 t^[0,T] 




= 0 . 


(3.28) 


(3.29) 


Proof. See Appendix [Ij 


Since the payoff of a typical FX or equity contract grows at most linearly in the exchange 
rate or the asset price, we only need to know the strong convergence in of the discounted 
process to deduce the convergence of the time-discretization error to zero. The following 
theorem can be generalized to the case relatively easily, for all cr > 1, upon noticing that 


the critical time T* from (3.16) - (3.19) is always greater than the one from (3.21) - (3.22). 


Theorem 3.10. If 2kf9f > and T < T*, then the diseounted proeess eonverges strongly 
in L^, 


i.e., 


where T* is given below: 


lim sup E 

St^o te[o,T] 


\Rt — Rt 


= 0 , 


T* = < Psv^ - k 


, if k < psvi 


(3.30) 


(3.31) 


oc , if k> psv^> 


Proof. The convergence in probability of the discounted process is a consequence of Propo¬ 
sition by taking the domestic interest rate to be zero. The rest of the proof follows the 
argument of Theorem |3.9| closely and makes use of Propositions |3.4| and |3.5[ ■ 


We can extend the convergence analysis from the four-dimensional Heston-CIR model to 
multi-factor short rates with CIR dynamics and a term structure, in which case Propositions 
3.4| to [3^ and hence Theorem |3.10| still hold, albeit with slightly modified proofs. 


The condition k > Psvi^ also known as the good eorrelation regime (Jacquier and Martini 


2011), is almost always satisfied in both FX and equity markets. This is because the speed of 


mean reversion k is usually larger than the volatility of volatility (^, a fact clearly illustrated 
in Table And even if this was not the case, a negative correlation between the underlying 
process and the variance, as is typically the case in equity markets (the so-called leverage 
effeet)^ or a small absolute value of this correlation, as is typically the case in FX markets, 
would ensure the validity of the condition. Furthermore, the Feller condition 2kf9f > in 
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Table 1: The calibrated Heston parameters. Column 2: for USD/EUR market data of 2 January 
2004 - 27 September 2005 ( Jessen and Pouls^|2013 ). Column 3: for EUR/USD market 
data of 22 August 2006 (Elices and Gimenez 2013[ ). Column 4: for the S&P 500 index 
between 2 January 1990 - 30 September 2003, using VIX data (Ai’t-Sahalia and Kimmel 
Column 5: 


2007[ ). Column 5: for the S&P 500 index between 2 January 1990 - 30 December 2011, 
using two out-of-the-money options written on the index ( Hum et al.||2014 ). 


Parameter 

Jessen & Poulsen 

Elices & Gimenez 

A’it-Sahalia & Kimmel 

Hum et al. 

k 

2.2200 

1.1000 

5.1300 

1.9775 

e 

0.0120 

0.0097 

0.0436 

0.0376 

e 

0.1830 

0.1400 

0.5200 

0.4568 

Psv 

0.0634 

0.1400 

-0.7540 

-0.7591 


( |3.28[ ) for the foreign interest rate is generally satisfied in practice, a fact clearly illustrated 
in Table 1^ We do not require a Feller condition for the stochastic volatility, and this is not 
always given in practice. 


Table 2: The calibrated Cox-Ingersoll-Ross parameters. Column 2: to the 3-month US Treasury 
bill yield between January 1964 - December 1998 ( Driffill et al.||2003 ). Column 3: to the 
US Treasury bill yield between October 1982 - April 2011 (Erismann 2011). Column 4: 
to the Euro ATM caps volatility curve on 17 January 2000 (Brigo and Mercurio 2006). 
Column 5: to the Euro OverNight Index Average between 1 January 2008 - 6 October 
2008 ( Laffers|[2009 ). Column 6: using historical data for Euro between 1 January 2001 - 
1 September 2011 ( Amin||2012 ). 


Parameter 

Driffill et al. 

Erismann 

Brigo & Mercurio 

Laffers 

Amin 

^dj 

0.0684 

0.1104 

0.3945 

0.2820 

0.1990 

^dj 

0.0161 

0.0509 

0.2713 

0.0411 

0.0497 

^dj 

0.0177 

0.0498 

0.0545 

0.0058 

0.0354 


3.3 Option pricing 

We conclude this section with a brief study on the convergence of mixed Monte Carlo/PDE 
estimators for computing EX option prices. Define the fair price of an option written on S: 


U = E 




and its approximation under (2.11): 


[/ = E 




(3.32) 


(3.33) 
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where the payoff function / may depend on the entire path of the process and all expecta¬ 
tions in this section are under the domestic risk-neutral measure Q. The following theorem 
is concerned with the convergence of the time-discretization error to zero under the four- 
factor FX model, and can be extended to multi-factor CIR short rates in a straightforward 
manner. The proof employs Propositions |3.4[ |3.5| and 3.8, and Theorem 3.10[ However, we 


omit it here because of its similarity to the proof of Theorem 2.1 in Cozma and Reisinger 


(2015a), once the aforementioned results are established. 

Theorem 3.11. Suppose that ‘2kf6f > Then the following two statements hold: 

(i) The approximations to the values of the European put, the up-and-out eall and any 
harrier put option defined in (3.33) eonverge as 5t ^ 0. 

(ii) IfT < T*, with T* from ( |3.31| ); the approximations to the values of the European eall, 
Asian options, the down-and-in/out and the up-and-in eall option defined in (3.33) 
eonverge as 5t ^ 

For European contracts, we can evaluate the conditional option price, i.e., the innermost 
expectation in (2.13), analytically. Hence, consider M simulations of the discrete paths of 
the variance and the interest rates and, for 1 < j < M, let ouj denote the j-th sample. Then 


1 

M 


M 


e-/o 


-JAv _ 


(jJ = (jJn 




= I af(3.34) 


M 


M 




is the mixed Monte Carlo/PDE estimator of the European option price at time t = 0. The 
global error can be split into two parts: 


Error = E 


M 


e-/o ^ E ^ I = cvj 

i=i 


= E 


e-So^tdiff^Sr) -E e-/o>t'^V('5T) 




+ E 


E 


e- fo ^tdtf^Sr) I | 


M 




M 




The first term is the time-diseretization error and the second term is the statistieal error. 


The convergence to zero of the former was derived in Theorem 3.11 for European put and 
call options, result that can be easily extended to other financial derivatives, such as binary 
options. The convergence to zero of the latter follows from the Central Limit Theorem (see 
Glasserman||2003 ) upon noticing the following upper bound on the variance. 


Var ( E 


o-fo 


/(5r)|a 


f,d,v 

T 


) < E 





O-fo 


f{ST)\g^, 


f,d,v 

T 


< E 


.-2/o rfdt 


fiSrf 
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Assuming / has at most polynomial growth, we can employ Proposition |3.3| to deduce the 
finiteness of the variance under some conditions on the model parameters. In particular, if 
/ is Lipschitz, Proposition |3.5| gives sharper sufficient conditions, for ol — 2. 

For path-dependent contracts, closed-form solutions are rarely available and we rely on 
finite differences instead to compute the conditional option price. Conditioning on the j-th 
realization of the variance and interest rates paths, let Si) and St\ P, L) be 

the solutions to the conditional PDF and to the associated finite difference scheme, when 
a uniform mesh with P time steps and L spatial steps is employed, respectively. It can be 


shown that the conditional PDE has a unique solution (see Section 7.1.2 in Evans 1998) 


which is, in fact, the conditional option price (see Theorem 7.3.1 in Shreve|[2004 ). Then 


M 

M E «“”(»■ So;-p.^) 


(3.35) 




is the mixed Monte Carlo/PDE estimator of the option price at t = 0. The global error can 
be split into three parts: 


Error = E 


M 


i=i 


= E 




-E 


+ E 


E 


e-fo^tdtf^S) I E ^ e-fo^dtf^s) \ = a;. 


e-/o 

1 ^ 


M 




J = 1 


1 ^ 


J = 1 

The first term is the time-discretization error^ the second term is the statistical error and the 
third term is the finite difference (FD) discretization error. The convergence to zero of the 
first term was derived in Theorem 3.11| for Asian and barrier options, and the convergence 
to zero of the second term is a consequence of the Central Limit Theorem (see the above 
discussion). However, the convergence of the third term, i.e., of the finite difference scheme, 
depends on the contract and the particular scheme employed. 


4 Variance reduction analysis 

In this section, we carry out a variance reduction analysis for European option valuation 
with the mixed Monte Carlo/PDE method under the four-factor EX model. However, the 
theory extends naturally to general interest rate dynamics. We use standard Monte Carlo 
with the log-Euler discretization as the reference method and define X and S to be the 


time-continuous approximations of x, defined in (3.24), and defined in (2.1), respectively. 
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Then 


Xt = X. + 


(< - < - 5 Vt.) (t - «„) + vt AIV,‘, 


where AW/ = W/ — W/^ whenever t € Integrating (4.1) leads to 


St = -So exp 


— ff — 

' u ' u 




du T 




(4.1) 


(4.2) 


We prefer the log-Euler scheme to the standard Euler scheme because it preserves positivity. 
Moreover, if the processes r^, are constant, then the first scheme is exact. Recall that 


St = 5oexp 




ft [ZZ ^ 

du + ail / Y 

Jo 


^ r 

/ 

i =2 ■^0 


E 72 ^ du 

ot 


Since E is piecewise constant, we deduce that St^ = 5^^, VO < n < A^. The quantity that 
we want to estimate is the fair price of a European option with payoff function /, i.e.. 


0 = E 


e-/o "*^V(Sr) 


(4.3) 


Then the corresponding standard and mixed Monte Carlo estimators are 

M 

M 


1 ^ T 

©stdMC ^ T7 E ^‘^■^''^V(ST(i)), 


i=i 

M 


OmixMC 


J = 1 


Define 


and 


Var,tdMC = Var(0,tdMc) = ^Var (e-^o""‘‘'V(ST)) 
Var^ixMC = Var(0™,Mc) = ^Var (e [e-/(S t) \ 


(4.4) 

(4.5) 

(4.6) 

(4.7) 


Let the variance reduction factor (as in |Dang et al.|[20T5 ) and the standard deviation ratio 

Vax^t^ and Tdev = (4-8) 

VarmixMC 

For convenience, we also define the discount factors, D = e~ Jo and D — e~ Jo 
Remark 4.1. From the “tower property” of conditional expectations, since St = 

IE [©stdMc] = E [DfiST)] = E [DfiSr)] = E [E [DfiSr) \ G^;^’^]] = E [0n,ixMc] • 
Therefore, the standard and the mixed Monte Carlo estimators have the same discretization 


17 















bias, i.e., 


(4.9) 


Bias(0stdMc) —Bias (©mixMc) • 

Remark 4.2. From the law of total variance we know that 

Var(5/(5r)) = Var(E [.D/(5r) | + E [Var(5/(5r) | . 


(4.10) 


However, since St = St and the variance is non-negative, we deduce from dlel) - dTH) that 
the variance of the standard Monte Carlo estimator is greater than or equal to the variance 
of the mixed estimator, i.e., 

Varg^jMQ ^ VarjnixMC • (^*11) 


Assuming a non-trivial payoff function /, equality occurs in (4.11) if and only if the second 
expectation on the right-hand side of (4.10) is zero, i.e., if and only if St is ^^’^’^-measurable. 
Since S = AA^ and after some straightforward calculations, we find an equivalent condition: 

an = 0 ^ l- pla- plf - plf + 2p^dpvfpdf = pIv (l “ plf) + pL(1 “ plf) + plf (l “ pld) 
4“ 2psvPsd(^Pvf Pdf Pvd^ 4“ 2psvPsf (^PvdPdf Pvf) 4“ 2pgdPsf (^PvdPvf Pdf)- 
In particular, if the variance and the two interest rates are pairwise independent, then 

an = 0 + plf = 1. 

Therefore, apart from this case, combining conditioning with Monte Carlo always reduces 
the variance of estimates. This is, however, to be expected since we eliminate the additional 
noise that comes from simulating the Brownian motion W^. 


Remark 4.3. Note that for any 2 < j < 4, using the Ito isometry, we have 

rT /— A VrT 

dt 


Var u 


(4.12) 


Moreover, using Cauchy’s and Holder’s inequalities, Fubini’s theorem, and Remark 3.2 and 
Propositions 3.4 and 3.5 in Cozma and Reisinger (2015a), one can easily prove that 


E 

lim - 

5t^o Var 




E 


!q vtdt 


(if Vtdt 


Var 


{lo Vtdt 


(4.13) 


From Dufresne (2001), we can compute the first two moments of the integrated square root 
process explicitly. Hence, we find 

rT 


E 


/ 

L^o 


vt dt 




(4.14) 
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and 


Vai-I 1^ v,dtj = p 




Vtdt 


<|,(l + e--)E[^ 

<^{l + e->‘^)Ey\,dt . 


+ «! 

fe 2 


Ae--(4 + 2 tT 


.-kT 




(4.15) 


Combining (4.12) - (4.15), we deduce that for sufficiently small values of 6t^ 


Var WtdWl 


Var 


(/o" 


Vfdt 


> 


> 


^ 2 (i + g-fcT) 2 ^ 


(4.16) 


The data in Tables and suggest that the interest rates have little impact on the variance 
of the mixed Monte Carlo estimator, and also that fc ^ in both FX and equity markets. 
Hence, the stochastic integral on the left-hand side of (4.16) - part of the dividend yield q 


defined in ( 2 . 12 ) - contributes to the overall variance of the mixed estimator considerably 
more than the squared volatility defined in (2.12). Therefore, we expect the minimum 


variance of the mixed estimator to be attained when all but the first term in q vanish, i.e., 
when ail = 1 and ai 2 = ai 3 = ai 4 = 0 . 


In practice, the volatility of interest rates in the CIR model is very small, i.e., ^ 1? 

a fact which can be observed in Table Moreover, the volatility of volatility in the Heston 
model - calibrated to FX market data - is also small, i.e., <C 1, and is significantly smaller 
than the rate of mean reversion, i.e., <C A:, a fact clearly illustrated in TableHence, the 
drift term in the square root model for the variance dominates the diffusion term. Therefore, 
we assume in the subsequent analysis “almost deterministic” dynamics of the variance and 
interest rates. Let 7 , 7 ^ and 7 -^ be the FTE discretizations of 'u, and corresponding 
to = 0, i.e., when the volatility of volatility parameters are equal to zero. Then 


[ rfdt^ [ ^fdt, [ r{dt^ [ ^{dt, [ Vtdt^ [ jtdt, [ JvtdW^^ [ ^tdWl 

JO JO JO JO JO JO JO JO 

Remark 4.4. Suppose that is independent of the Brownian motions VF^, and VF-^, 
i.e., that an = 1 and ai 2 = ai 3 = ai 4 = 0. Employing (4.7) as well as the above assumption 
on the dynamics of the variance and the domestic and foreign interest rates, 

VawMC = ^Var (e”E [/{St) \ ~ ^e-2/o^^*‘^War(E [/(5r)]) = 0. 
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On the other hand, since St = St^ we know from (4.6) that 


VarstdMC = (e -^o Var(/(S'T))• 

Assuming a non-trivial payoff / as well as a non-zero squared volatility v, since the variance 
of the mixed estimator is close to zero, this results in a substantial variance reduction. 


Remark 4.5. Let / be the European call option payoff and suppose that the dynamics of 
the variance and interest rates are “almost deterministic” (i.e., their own volatility param¬ 
eters are small compared to their mean-reversion speed). Then we can approximate the 
discounted payoff as follows: 



For now, assume that an > 0. For convenience, define the quantities 


(4.17) 


(4.18) 


as well as 




and b = 


log{F/K) + (0.5 + ^2)^2 




(4.19) 


a 


We use (4.17), the conditional option price formula in (2.14) and differentiate the variance 
of the mixed Monte Carlo estimator with respect to the parameter g defined in (4.18) to 
find, using lengthy integration by parts. 


dg 


Var^ixMC - 


A E [$(aZ + 6)2], 


(4.20) 


where $ is the standard normal CDF and Z ~ V(0,1). Note that if an < 1, i.e., if ^ > 0, 
the right-hand side of (4.20) is strictly positive. This implies that the variance of the mixed 
estimator decreases as an increases, and attains its minimum when an = 1. In fact, we 
know from Remark 14.41 that 

VarmixMc(«ii = 1) « 0. (4.21) 


Moreover, 


E [4>(aZ + 6)2] = $ 


vT+~c? 


2T 


1 


\ vT+~c? 


Vl + 2a2 


(4.22) 
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where Owen’s T function ( Owen|[r980| ) is 

TiP,!)) = (p{(3) f with T(^,l) = _ 1$(/5)2^ 

Jq 1 X z z 

and 6 is the standard normal PDF. Note that 


(4.23) 




1 1 F 2^~ 

— — log ^ o ^ jO' and 


a 


K 


' 1 — 


Vl + 2a2 


=01 


=02{S) 


' 1 + ^' 
=0is) 


,2 ■ 


(4.24) 


According to the data in Table [T] 0^1 and so 7 ^ 1 . Hence, for low maturities T we have 
d ^ 1. From Remark 4.2, using ( |4.21 ), we know that the variance of the standard Monte 
Carlo estimator is obtained by integrating (|4.20|) over [0,1]. Therefore, 


^ /oV[^(/3i + /32(0) -2T(/3i+/32(0,^^(0)]dt/ 

So 'A.^{h ++ P2{i'),^{i'))]du' 

Ffom ( |4.23[ ) and ( |4.24[ ), we deduce that € [0,1], 

$(/3i + 0.5d)^ < $(/3i + /32(0) “ 2T(/3 i + P2{i^)','d{y)) < ^(/3i + 1.5d). 

Therefore, 


4>(/3i + 0.5d)' 


1 / -n / ^(A + l-5d) 

2 — ^ var ^ 


(4.25) 


(4.26) 


(4.27) 


4>(/3i + 1.5d) $(,5i + 0.5d)" Q" 

The inequalities in (|4.27|) are approximate in the sense that the variance reduction factor is 


bounded from above and below by approximated quantities. Assuming that d <C 1, (4.27) 
becomes 


^(/5i) 


^ hyar ^ 


1 


$(/3l)p2 


(4.28) 


Moreover, if we also assume that the option is deep in-the-money, then /3i is relatively large 
(e.g., /3i > 1 ) and so y^$(/3i) 1 . Hence, we conclude that 


Tde 


(l - aji) 2, Van > 0. 


(4.29) 


However, we know from Remark 4.2 that the standard deviation ratio is one when an = 0. 


Hence, (4.27) - (4.29) hold for all values of an. Interestingly, the deep in-the-money case 


is also the one where the payoff can be treated as smooth, and this also explains Figure]^ 
Note that Remark 14.51 is consistent with Remarks 14.21 
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5 Numerical results 


In this section, we test the efficiency of the mixed Monte Carlo/PDE method by comparison 
with alternative numerical schemes for two derivatives: a European call and an up-and-out 


put option. For the former, we use an analytical formula (2.14) to compute the conditional 


option price and benchmark against standard Monte Carlo with a log-Euler discretization 
and, under a simple correlation structure, against the semi-analytical formula of |Ahlip and 


Rutkowski (2013). For the latter, we use either finite differences or the perturbative formula 


of Fatone et al. (2008) for the conditional option price, and Monte Carlo - with or without 


Brownian bridge - as the reference method. 

First, we demonstrate the convergence of the mixed Monte Carlo/PDE method and find 
empirical convergence rates. Then, we investigate the sensitivity of the variance reduction 
factor ( |4.8[ ) to changes in the model parameters and link the numerical results to the analysis 
in Section Our platform for the numerical implementation was MATLAB 2012a. The 
machine configuration on which all numerical tests were conducted is: Intel(R) Core(TM) 
i3 CPU, M370, 2.40 GHz, Memory (RAM): 8.00 GB, 64-bit Operating System running 
Windows 7 Professional. 

Throughout this section, we assign the following values to the underlying model param¬ 
eters as a base case, and vary a selection individually or jointly: = 0.0275, Tq = 0.0524, 

rl = 0.0291, k = 1.70, kd = 0.20, kj = 0.32, 6 = 0.0232, 0d = 0.0475, 0f = 0.0248, 
^ = 0.1500, = 0.0352, = 0.0317, = -0.10, p,d = -0.15, p^f = -0.15, pw = 0.12, 

Pyf = 0.05 and pdf — 0.25. These values are consistent with empirical observations in FX 
markets and are close to the calibrated values in Tables 0 and Also, the values of the 
correlations psf and pdf are borrowed from Piterbar^ ( |2006[ ). 


5.1 European call option 

Let the spot, the strike and the maturity be: Sq = 105, K = 100 and T = 1.5. On a side 
note, FX option quotes are in terms of volatilities for a fixed delta and a fixed time to expiry, 
and not in terms of strikes. However, the strike price corresponding to a quoted volatility 
can easily be recovered using a conversion formula ( Ahlip and Rntkowski|2013 ). Recall that 


0 from (4.3) is the true option price, whereas ©stdMC from (4.4) and ©mixMC from (4.5) are 


the standard and the mixed Monte Carlo estimators, respectively. In Table we report the 
common discretization bias, i.e., the time-discretization error, which is the same for both 


estimators by Remark [4.If and the two standard errors, i.e., the two standard deviations of 
the sample means, which we denote by StDev and estimate using 10000 samples. However, 


according to Ahlip and Rutkowski (2013), a closed-form solution for the European option 


price is not available under a full correlation structure. Hence, we need to find an accurate 
reference estimate ©* in order to evaluate the bias. We employ the mixed algorithm with 
M — 2x10^ simulations and N = 200 time steps to find: ©* = 12.11968. Then, the bias for 
a specific number of time steps is estimated using the reference estimate ©*, whose accuracy 
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is discussed below, and a sufficiently large number of simulations, i.e., M — 2x10^. 

The time needed to obtain a call price estimate with 64000 simulations and 8 time steps 
is 0.25 seconds for the standard Monte Carlo method and 0.22 seconds for the mixed Monte 
Carlo method. Hence, the computational cost is 12% lower with the latter. This, however, 
is to be expected since we only simulate three of the four underlying processes. 

Table 3: Simulation results for a European call option. 


Time steps 

Discretization bias 

Simulations 

StDev (stdMC) 

StDev (mixMC) 

1 

0.37231 

1000 

0.50348 

0.09342 

2 

0.08039 

4000 

0.24126 

0.04194 

4 

0.01775 

16000 

0.12042 

0.02023 

8 

0.00444 

64000 

0.06071 

0.00994 

16 

0.00160 

256000 

0.02932 

0.00487 

32 

0.00073 

1024000 

0.01507 

0.00262 


The data in Table indicate that the bias is small even when only a few time steps are 
used. This phenomenon can be explained by the fact that the Feller condition is satisfied, 
i.e., 2k9 > 2kd0d > ^ C/? hence the discretizations of the variance and 

interest rate processes rarely hit zero. The simulation results in Table confirm the square- 
root convergence of the statistical error and the first-order convergence of the discretization 
error. Then, using extrapolation, we obtain an approximate root mean square error (RMSE) 
of the reference estimate: 


Bias(0*) « 1.168x10“^, StDev (0“^) « 0.593x10 


-4 


RMSE(©*) 1.310x10 


-4 


This is equivalent to an RMSE of about 0.001% of the actual option price, suggesting that 
the reference estimate 0* = 12.11968 is accurate to three decimal places. 

Furthermore, the standard deviation is reduced by 83%, i.e., is approximately six times 
lower with the mixed method. However, the variance reduction is not consistent across the 
range of possible values of the model parameters, and is most sensitive to changes in the 
correlations between the exchange rate and the squared volatility or the interest rates, i.e., 
Psv^ Psd si^nd psf^ the speed of mean reversion k and the volatility of volatility 

The values pyd — 0.12, p^f = 0.05, pdf = 0.25 and psd^ Psf ^ [—0.5,0.5] lead to 
a valid (i.e., positive definite) correlation matrix. The data in Figure suggest that the 
highest variance reduction is achieved when the absolute values of the correlations psd 
and Psf are small or, more precisely, when psv ~ —0.05 and psd ~ Psf ~ 0, in which case the 
standard deviation is reduced by a factor of 20. Hence, the best performance of the mixed 
estimator coincides with the independence of from the Brownian motions and 


wf, an observation which is consistent with Remarks |4.3| and |4.4| 


Furthermore, we infer from Figure that, in addition to the lower computational cost, 
the mixed Monte Carlo/PDE method outperforms standard Monte Carlo in terms of accu- 
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Figure 1: The standard deviation ratio Tdev from (4.8) with 4000 simulations and 10 time steps, 
plotted against the correlation coefficients psd and psf when the last two are equal. 


racy for all psd-) Psf ^ [—0.5,0.5], and that the variance reduction factor decreases as 
the absolute values of the correlations between the exchange rate and the squared volatility 
or the interest rates increase, two facts that were recognized in Remark |4.2[ 

Assuming that psd — Psfi since p^d ~ Pvf ^ 0, we may approximate 

piv-lp'id- ( 5 - 1 ) 

However, 

Tdev « (1 - di)"" ~ (pL + 1-6pL)”" (5.2) 

from Remark 14.51 Combined with the above observation on the location of the maximum 
variance reduction factor, we infer that the set of points {psv^ Psd) corresponding to T^ev = 
for some /i > 1, takes approximately the form of an ellipse with the equation: 

{p,, + 0.05)^ + = 1- (5-3) 

This confirms that the isolines in Figure [^display an elliptic shape, a fact also illustrated in 
Table 1^ where the correlation pairs (psv^Psd) are chosen so that p — 2.812. The estimated 
standard deviation ratio along the contour line is 3.148 ± 0.239, which is close to the theo¬ 
retical value p. This attests to the accuracy of our approximations and numerical results. 

Next, we assume a partial correlation structure between the Brownian drivers such that 
the squared volatility dynamics are independent of the domestic and foreign interest rate 
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Table 4: The estimated standard deviation ratio for different correlations. 


Psv 

-0.40 

-0.40 

-0.30 

-0.30 

0.20 

0.20 

0.30 

0.30 

Psd 

-0.05 

0.05 

-0.20 

0.20 

-0.20 

0.20 

-0.05 

0.05 

Fdev 

3.288 

3.088 

3.387 

2.910 

2.951 

3.127 

2.952 

3.055 


dynamics, i.e., = p^f — 0. Using Remark 4.5 and computing an explicitly from (2.2), 

1 


Tdev 


pL + 


1 - 


Pdf 


{Psd + Psf ‘^PsdPsfPdf) 


(5.4) 


Before, we fixed p^f and analyzed the variance reduction with respect to psv, Psd ^'iid Psf, 
when the last two were equal. Now, we fix psv instead and focus on the effect of varying pdf 
on Tdev To this end, one can easily show that 


Psd T Ps/ ‘^PsdPsfPdf ( 2 2 \ 

-—^- > max [p^d^ p^p. 

^ Pdf 


(5.5) 


Assuming that psd and psf are not simultaneously zero, equality in (5.5) holds when 


_ * _ Psf - Psd 


\Psd\^\Psf\ 


(5.6) 


Upon its substitution into (5.4), we find a theoretical standard deviation ratio p. Table 
estimates T^ev and compares it with /i, when pgy — —0.10, Psf = —0.20, Psd ^ [—0.25,0.25] 
and Pdf G [—0.85,0.85]. Note that these correlation values guarantee a symmetric positive 
definite correlation matrix. Hence, when lies outside the interval of allowed values, we 
estimate T^ev using pdf = ±0.85 instead. The data in Table suggest that the approxima¬ 
tion (5.4) to the standard deviation ratio is quite accurate, especially for a strongly negative 
correlation between the two interest rates. 

Table 5: The empirical and theoretical standard deviation ratios for different correlations. 


Psd 

-0.25 

-0.20 

-0.15 

-0.10 

-0.05 

0 

0.05 

0.10 

0.15 

0.20 

0.25 

Pdf 

0.80 

1.00 

0.75 

0.50 

0.25 

0 

-0.25 

-0.50 

-0.75 

-1.00 

-0.80 

Fdev 

4.501 

5.301 

5.305 

5.224 

5.073 

5.010 

4.850 

4.867 

4.554 

4.467 

3.788 

P 

3.714 

4.472 

4.472 

4.472 

4.472 

4.472 

4.472 

4.472 

4.472 

4.472 

3.714 


By close inspection of the data in Figure]^ we infer that for each psd G [—0.25,0.25], the 
highest variance reduction is achieved for as defined in (5.6). Hence, we conclude that 
p exhibits the qualitative behaviour of F^ev? so that (5.4) provides a good approximation 
to the standard deviation ratio. In fact, our observations suggest that p acts as a lower 
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bound. We can extend these results to a full correlation structure between the Brownian 
drivers as long as and pyf are close to zero, as seen before. 



Correlation (s,d) 


Figure 2: The standard deviation ratio with 1000 simulations and 4 time steps, plotted against the 
correlation coefficients psd and pdf when psv = —0.10, psf = —0.20 and p^d = Pvf = 0- 


Suppose that the exchange rate dynamics are independent of the interest rate dynamics, 
i.e., that psd — Psf — 0, and define the optimal correlation to be the value corresponding 
to the maximum standard deviation ratio, for a given volatility of volatility Then the 
data in Figure [^suggest that the highest variance reduction is achieved when the absolute 
value of the correlation is small. In fact, we notice two things. First of all, that the optimal 
correlation approaches zero as the volatility of volatility decreases, i.e., that lim^^o P'lv — 0? 
which is to be expected from Remark |4.3[ And second, that the standard deviation ratio 
at increases as the volatility of volatility decreases. In practice, ^dj 1 (see Table [^, 
so the two interest rates have little impact on the variance of the mixed Monte Carlo/PDE 
estimator. Therefore, since is close to zero, the variance comes mainly from a defined in 
(2.12). On the other hand, smaller values of ^ result in a smaller variance of a, and hence 
of the mixed estimator as well, which leads to a higher standard deviation ratio. 

We observe a similar behaviour when increasing the speed of mean reversion k or the 
long-run variance 9. Larger values of k result in a smaller variance of a because of the mean- 
reverting property of the squared volatility, which ensures that the process returns to the 
long-run average quickly. On the other hand, larger values of 9 produce higher volatilities, 
which then leads to an increase in the variance of the standard Monte Carlo estimator due 
to the larger diffusion term in the SDE driving the exchange rate process. We conclude the 
analysis by noting that values of psd and psf close to zero produce similar results, to some 
extent. However, when the absolute values of the two correlations are not small, the impact 
of (7, and hence of fc, ^ and 0, on the variance of the mixed estimator is reduced. 

The maximum in Figureis attained around psv = 0 and ^ = 0.05, where the standard 
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Figure 3: 


The standard deviation ratio with 4000 simulations and 10 time steps, plotted against 
the correlation psv and the volatility of volatility ^ when psd = psf = 0- 


deviation ratio is Tdev = 23. Therefore, the same level of accuracy with the standard Monte 
Carlo method requires 529 times more simulations. 

Next, we examine the variance reduction for different spots and maturities. The data in 
Figure suggest that unless the option is far out-of-the-money and the maturity is small, 
varying *So and T has little impact on the standard deviation ratio, which is approximately 
2. Computing the option price with standard Monte Carlo means integrating the payoff 
function, which is not differentiable at the strike, whereas with the mixed Monte Carlo/PDE 
method we integrate the smooth conditional price. As the probability of a positive payoff 
decreases with Sq and T, estimating it accurately with the former requires more simulations. 
Hence, the relative standard error of the standard Monte Carlo method increases as we go 
farther out-of-the-money or approach maturity, and the benefit of employing the mixed 
algorithm becomes clear. For example, for a 3-month call with a spot at 75% of the strike, 
we observe a variance reduction factor of 5275. 

Furthermore, Figure [^suggests that the approximation for the standard deviation ratio 
derived in ( 4.29[ ) does not hold when the option is far out-of-the-money. For instance, we 
compute ail — 0.7893, which gives a theoretical standard deviation ratio pi = 1.629. When 
Sq — 105 and T = 0.5, this is close to the estimated value F^ev = 1.832. However, when 
Sq = 85 and T = 0.1, we observe a much higher standard deviation ratio F^ev = 19.863. 

Suppose that the domestic and the foreign short rate dynamics are independent of each 
other, and also independent of the dynamics of the exchange rate and its volatility, and let 
Psv = —0.10, as before. Moreover, suppose that the other model parameters, as well as the 
spot, the strike and the maturity, take the values listed at the beginning of this section. We 
use the mixed Monte Carlo/PDE method with M = 8x10^ simulations and N = 200 time 
steps, for an RMSE of about 0.001%, to obtain a call price estimate: 0' = 12.13621, which 
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Figure 4: The standard deviation ratio with 1000 simulations and 4 time steps, plotted against the 
spot So and the maturity T when psv = —0.60. 


is about 0.14% higher than the estimate corresponding to a full correlation of the factors, 
i.e., 0* = 12.11968. On the one hand, the postulated independence of the factors is critical 
from the point of view of analytical tractability ( Ahlip and Rutkowski||2013 ), but can result 
in fairly different option prices. On the other hand, a full correlation structure leads to a 
richer model and a better fit to the observed market data. Finally, we test the accuracy of 
the mixed method and employ the semi-analytical pricing formula of | Ahlip and Rutkowski 


(2013) to find the true option price, 0 = 12.13603, and thus a relative error of O' of about 


0.0015%, which confirms that mixed Monte Carlo/PDE estimates are correct. 

5.2 Up-and-out put option 

Let the spot, the strike, the barrier and the maturity be: So — 100, K = 105, B = 110 and 
T = 0.25, and consider a continuously monitored up-and-out put option. We will first value 
the contract using the mixed Monte Carlo/PDE method. Hence, for a specific realization 
of the variance and interest rates paths, we compute the conditional option price, i.e., 


u[ 


,(t.x) = E (A- - Sr f I 5“'. S, = x 

We know that u satisfies the following initial boundary value problem: 

dtu -h fitxdxU -h ^(i\iVtx^dxx'^ — = 0, \/0 < X < t < T 

u{t, B) = 0, yt < T 
u{T,x) = {K- x)+, \/0<x<B. 


(5.7) 


(5.8) 
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We solve the PDE backwards in time from the initial condition, on a rectangular domain 
with t G [0, T] and X G [O.TSq^B] that is discretized on a uniform grid with + 1 temporal 
nodes and L + 1 spatial nodes. We selected this particular lower boundary of the spatial 
computational domain to reduce the number of spatial nodes, at the same time making sure 
that the truncation error arising from our choice of the domain is negligible. We employ a 
central difference scheme to approximate the spatial derivatives and a linearity boundary 
condition ( Tavella and Randall|2000 ) stating that dxxU = 0 at the lower boundary, where the 
option is deep-in-the-money and the price can be regarded as linear in x. For convenience, 
we have employed the same time grid used for the discretization of the squared volatility and 
the interest rates. The final estimate of the barrier option price is a Monte Carlo average 
over a sufficiently large number of discrete trajectories of the Brownian motions 
and W^. 

Alternatively, we can use the perturbative formula of Fatone et al. (2008) to approximate 
the conditional option price, and then a simple Monte Carlo average to estimate the outer 
expectation. Hence, we call this numerical scheme the mixed Monte Carlo/Pert method. 


Fatone et al. (2008) approximate the up-and-out put option price in the Black-Scholes 


model with time-dependent parameters via a series expansion and provide explicit formulae 
for the first three terms, which involve some elementary and nonelementary transcendental 
functions. However, we will focus only on the zeroth-order approximation because using a 
first-order correction term results in a hundredfold increase in computation time, and hence 
in a poor performance of the scheme as opposed to the mixed Monte Carlo/PDE method. 
Using M = 4x10^ simulations and N = 200 time steps to minimize the sampling and the 
discretization errors, respectively, we obtain a zeroth-order approximation: = 5.7700. 

Since a closed-form solution to the option pricing problem is not available, we need to 
find an accurate reference estimate ©* for the true option price © in order to compute the 
different errors of the numerical methods. Therefore, we use the mixed Monte Carlo/PDE 
algorithm with the Crank-Nicolson scheme, with M = 4x10^ simulations, N = 200 time 
steps and L = 20 spatial steps, to find: ©* = 5.7631. Hence, the approximation error of 
the mixed Monte Carlo/Pert method is: ©^ — ©* = 0.0069, i.e., about 0.1%. 

In Figure we report the time-discretization errors computed using the reference esti¬ 
mate ©* - whose accuracy is discussed below - or ©^, and a sufficiently large number of 
simulations and spatial steps, i.e., M — 4x10^ and L — 20. For the standard Monte Carlo 
and the mixed Monte Carlo/Pert methods, the time-discretization error is defined as the 
bias, whereas for the mixed Monte Carlo/PDE algorithm, due to our choice of the finite 
difference grid, it contains the finite difference (FD) time-discretization error. Henceforth, 
the term “discretization error” stands for the time-discretization error. 

On a side note, crossings of the barrier are monitored only at discrete times by standard 
Monte Carlo. This gives rise to a monitoring error, which is included in the discretization 
error. Moreover, due to the knock-out feature of the option, the true price is smaller than 
the Monte Carlo estimate, which explains the strong positive bias displayed in Figure 

The data in Figure [^suggest a square-root convergence of standard Monte Carlo and a 
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Figure 5: The log-plot of the time-discretization errors of the mixed Monte Carlo/PDE, the mixed 
Monte Carlo/Pert and the standard Monte Carlo methods, against the time step. 


first-order convergence of the mixed algorithms. On the other hand, we can use a Brownian 
bridge technique (see Glasserman||2003] ) to improve the first method and recover the first- 
order convergence. Indeed, the (red) Monte Carlo curve in Figure]^ would almost coincide 
with the (green) mixed Monte Carlo/Pert curve with the Brownian bridge correction. For 
instance, we calculated the discretization bias with N = 8 time steps to be 0.0075 for both 
Monte Carlo with Brownian bridge and mixed Monte Carlo/Pert. 

In Figurewe report the space-discretization error computed using ©* and a sufficiently 
large number of simulations and time steps, i.e., M = 4x10^ and N = 200. The data suggest 
a second-order convergence as well as a local minimum discretization error when the strike 


price lies halfway between two adjacent nodes, a technique known as grid-shifting (Tavella 


and Randall 2000). Hence, considering the first-order convergence of the time-discretization 


error (T-Err) and the second-order convergence of the finite difference space-discretization 
error (S-Err) with the mixed Monte Carlo/PDE method, and using extrapolation, we obtain 
an approximate root mean square error (RMSE) of the reference estimate: 


T-Err 5.76x10 S-Err 10.80x10 StDev ~ 2.09x10 


-4 


1-4 


RMSE 1.67x10 


-3 


This is equivalent to an RMSE of about 0.03% of the actual option price, suggesting that 
the reference estimate 0* = 5.7631 is accurate to two decimal places. Next, we compare 
the three numerical methods in terms of computation time for a given level of accuracy, in 
particular, when the RMSE is at most 0.30% of the option price. First, using the empirical 
convergence rates determined above and extrapolation, we need M = 2.5x10^ simulations 
and N = 800 time steps, and hence a CPU time of 61.2 seconds, with the standard Monte 
Carlo method. Second, we reach this level of accuracy with the mixed Monte Carlo/PDE 
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Figure 6: The log-plot of the absolute space-discretization error of the mixed Monte Carlo/PDE 
method, against the spatial step. 


method when M = 12000, = 10 and L = 12, in 2 seconds. Third, we need to employ the 

mixed Monte Carlo/Pert method with a zeroth-order approximation, with M = 12000 and 
= 10, which takes 3.1 seconds. 

Therefore, when the np-and-ont put option price estimates need not be too accurate, e.g., 
when one decimal place of accuracy is sufficient, the two mixed algorithms are comparable 
in terms of CPU time and efficiency, and are considerably faster than standard Monte Carlo. 
However, a higher accuracy would require at least a first-order correction term in the mixed 
Monte Carlo/Pert approximation, making it highly time-consuming. We thus conclude that 
the mixed Monte Carlo/PDE method is the better of the three schemes. 

We mentioned above that Monte Carlo with Brownian bridge recovers the observed 
first-order convergence of the discretization error and the level of the bias from the mixed 
Monte Carlo/Pert method. The time-discretization error with the mixed Monte Carlo/PDE 
method is approximately 1.6 times larger, and includes the ED time-discretization error. 

For a two-decimal-place accuracy, we fix 100 time steps and 20 spatial steps, such that 
the space and time-discretization errors are about 0.02%. The time required to obtain a 
barrier option price estimate with 40000 simulations is then 26 seconds for the mixed Monte 
Carlo/PDE method and 2.1 seconds for Monte Carlo with Brownian bridge (1.4 seconds 
for standard Monte Carlo). Therefore, the computational cost is 92% lower with the latter. 
In conclusion, due to the square-root convergence of the standard deviation, T^ev needs to 
be above 4.5 in order for the mixed Monte Carlo/PDE method to outperform Monte Carlo 
with Brownian bridge. Just as in the European call option case, the variance reduction is 
most sensitive to changes in the correlations between the exchange rate and the squared 
volatility or the interest rates. 

Figure [^exhibits similar characteristics to those of Figure In particular, the highest 
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Correlation (s,v) 


Figure 7: The standard deviation ratio with 400 simulations, 10 time steps and 4 spatial steps, 
plotted against the correlation coefficients psd and psf when the last two are equal. 


variance reduction is achieved when psv ~ 0.05 and psd ^ Psf ~ 0, in which case T^ev = 40. 
Based on a previous observation, this is equivalent to a 80 times lower computational effort 
with the mixed algorithm. A careful inspection of the data in Figure suggests that the 
set of points (psv^Psd) corresponding to F^ev > 4.5 can approximately be described by the 
following inequality: 

(p,, - 0.05)'+ 1.6p2^< 4.2-2, (5 9^ 

i.e., the inside of an ellipse, this being the set of parameters where the benefit of variance 
reduction outweighs the additional complexity of solving the conditional PDF numerically 
in this instance. 

The variance reduction achieved by the mixed method results in computational savings 
(in the number of samples) by a factor roughly independent of the desired accuracy. Con¬ 
versely, higher accuracy of the finite difference method can only be achieved by a larger 
number of mesh points. Hence, it would appear that for high enough accuracy the mixed 
method can never win over the standard Monte Carlo method. An asymptotic complexity 
gain of the mixed method for small errors would require that the PDE can be solved with 
constant effort independent of the desired accuracy (and for this effort to be outweighed by 
the reduced number of samples required for a given statistical error). Multilevel Monte Carlo 
methods (Giles and Reisinger 2012) recently developed for stochastic PDEs are designed 
precisely for this goal, and achieve it by concentrating the dominant number of samples on 
the coarsest meshes, while computing corrections on finer meshes with a vanishing number 
of paths. The application and numerical analysis in the present context is the subject of 
further research. 
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6 Conclusions 


The numerical experiments carried out in Section suggest that the mixed method outper- 
forms both standard Monte Carlo and finite difference methods under certain circumstances 
depending on the contract and the model parameters. When a closed-form solution for the 
conditional option price is available, we usually see a considerable improvement in terms 
of both accuracy and computation time. If not, the mixed algorithm provides better accu¬ 
racy at the expense of added computation time and the set of parameter values where it 
outperforms the classical schemes is limited. 

The analysis carried out in this paper is not restricted to the four-dimensional Heston- 
CIR model, but can be extended to higher-dimensional problems. For instance, one may 
consider multi-factor short rates as in Dang et al. (2015), with CIR dynamics and a term 
structure, in which case the convergence and variance reduction analysis applies with some 
slight modifications of the proofs. Stochastic volatility accounts for volatility clustering, 
dependence in the increments and long-term smiles and skews, but gives rise to unrealistic 
short-term patterns in the implied volatility. Hence, in order to improve the behaviour of the 
implied volatility for short maturities, one could extend the original model to a stochastic- 
local volatility model as in Cozma and Reisinger (2015a), which can easily be implemented 
for barrier option pricing with no extra computational effort, or add an independent jump 
component to the spot FX rate, in which case analytical formulae may be available for the 
conditional prices of European options. For instance, that is the case when the distribution 
of the jump size is normal ( Merton||1976 ) or double-exponential ( Kon||2002 ). 

However, several unsettled questions remain, like the strong convergence rate of the dis¬ 
cretization scheme, or a finite difference scheme with an observed second-order convergence 
in time for pricing the barrier option. In addition, examining the hedging parameters is 
also relevant, and we intend to pursue all these topics in our future research. 
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Appendix A. Proof of Lemma 


3.1 


Let , 0 <t <T} be the natural filtration generated by and employ the shorthand 
notation [ • ] = E [ • . If we assume that t G [tn^ ^n+i] and condition on we get 


[®t] = 



Yu T 




exp 


A + 


t tji 2 



Upon noticing the identity below, 


sup Xx+l-iTx'^ = A1a>o, 

tcG[0,l] ^ 


we deduce that A < 0 implies E^^ [0^] < E^^ [©tn] A > 0 implies E^^ [0^] < E^^ [0t^+i]. 
Moreover, since Y is piecewise constant. 


/ 


in /^SWf 
6t 


du — 


ptn 

Jo 


Yudwy, V0<n<iV. 


Henceforth, we follow the argument of Proposition 3.6 in Cozma and Reisinger (2015b). 
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Appendix B. Proof of Proposition 


3.2 


We find it convenient to define a new stochastic process L by 
Lt = St exp 


J = Soexp |y (rt - ^Vu^du + J 


(B.l) 


Since St < Lt, Vt € [0,r], it suffices to prove the finiteness of the supremum over t of 

4 


E [L^] =S^K 


exp rtdu-^ [ Vudu + uj E J (®-2) 


Let 0 < t < T| be the natural filtration generated by the Brownian drivers and 

i.e., generated by the processes and v as observed until time T, and 0 < t < T} 
be the filtration generated by Conditioning the expectation on the right-hand side of 
(B.2) on the cr-algebra and taking into account that and are independent, we 
can compute the inner expectation using moment generating functions (MGFs): 


E 


exp Icjaii J vVI dW^ + ujai 2 J vVI dW^ 11 gf 
= E exp J vVI dw^ 1 1 E |exp I wai2 J s/v^dW^ 1 1 gf’ 

= exp { Y (oii + 012) y 


(B.3) 


Substituting back into (B.2) with (B.3) leads to 
E [L“] =S^E 


lcj‘^ (jjI ^ ^ pt .V 

exp<{w/ rf^du+ Y(“t + «i 2 )-x / Vudu+uj'^aij ^/n^ldWi'> . 

^ 0 _ _'^0 j 3*^^ 


Next, we employ Holder’s inequality with the pair (p, g), where p, g > 1 and ^ = p/(p — 1), 
in order to force the term J^r^du outside the expectation, and then condition the second 
expectation on the cr-algebra to arrive at 


exp {pcj J du 


E [L^] < 5 ^E 

-hgcc;ai4 / ^J^dW‘^ 
Jo 


E 


exp<g 


Y(a?i + a?2 + <Za?3)-| 


Vu du 


(B.4) 


All that is left to do is to show that the supremum over t of each of the two expectations 


on the right-hand side of (B.4) is finite. However, Proposition 3.2 in Cozma and Reisinger 
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(2015b) provides the following sufficient conditions: 

kd > 


as well as 


and 


k > qujpsvC 


+ \^<^Psv^k + [ah + ah) - q - k'^ <0, 


(B.5) 

(B.6) 

(B.7) 


for all T > 0. The first assumption in ( |3.7[ ) ensures that qo{a) > 1, and hence that qi{a) > 1. 
This implies that («)/(<?!(a) — 1) > 1. Due to the second assumption in (3.7), we can find 
p > 1 so that 

aqi(a) 


kj 

^>ap> —— -- 


kd > V^id- 


(B.8) 


The quadratic equation in x below has roots of different signs, with positive root go(<^)- 

+ \^oi,Psvik + (a\i + a\2) — — 0. 

However, the Holder pair satisfies p = q/{q — 1), so q < qi{a) < qo{a). Therefore, q lies in 
between the two roots of the quadratic, which implies that 


+ [‘^aPsv^k + [ah + ah) - - k‘^ < 0 . 


PYom (3.6), if psv > 0, 


q < qi[a) < 


(^Psv^ 


k > qapsv^, 


(B.9) 


(B.IO) 


and this clearly holds when the correlation coefficient is non-positive. Consider the three 
continuous maps below, which are strictly positive when oj = 

u^k- qupsvi\ uj^kd - 

u ^ k^ — u^h^ahq^ — \^a)psv^k + [ah + 012 ) “ 

Then we can find ai > a such that all three functions are positive on [a, ai). We have thus 
proved that conditions (B.5) - (B.7) are satisfied and the conclusion follows. The extension 
to the interval [l,(ai) follows immediately from Jensen’s inequality. ■ 

In the special case that ais = 0, i.e., psd = PsvPvd^ fhe argument is the same and the 
only difference appears in condition ( B.7[ ), which becomes 


[2upsv^k + (1 - ph) - k^ < 0 . 
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Henceforth, one can easily show that Proposition |3.2| still holds in this case as long as 


Appendix C. Proof of Proposition 


3.3 


For convenience, define a new stochastic process L by 

Lt = Stexpi^J = Soexp I du + an J \IVudWl 

* f77 5Wl 


j=2 *^0 


(C.l) 


As St < Lt, yt € [0,T], it suffices to prove the hniteness of the supremum over t and 5t of 
E [L“] = 


exp w / rfi^du—^^ Vudu + ojan / yVudW^ 
Jo 2 Jq Jq 


4 r 

+ ui'^aij / 
i=2 ■^0 


P?/, r. dU 
ot 


(C.2) 


Conditioning the expectation on the right-hand side on and bearing in mind that W^A. 
we can split the inner expectation into two parts, which we compute using MGFs. First, 


E 


exp (jjaii J yVu dW^ 


'7(1,V 


UJ 


= exp i f • 


f 


(C.3) 


Second, let t G [tn,tn+i)- As V is piecewise constant and has independent increments, 


E 


'■f 


. , 1= ml ^ 

exp <; ijjax 2 I Y 


7d,v 


= E 


rtn j y-_+ r^n+i j 

exp <j ojan ^ V ^ ‘ 


I 2 C t 7 J \ / 2 (^ - 

exp Y ^12 J ^udu I exp I — a 


/ 

Jtn 




r>d,v 

yT 


Vr,,du 


<exp{Y«i 2 y 


(C.4) 


Substituting back into (C.2) with (C.3) and (C.4) leads to an upper bound, 

4 


IE [A“] < Aq E exp |a; J r^du + J 


* 177 5Wl 


V,. du 
ot 
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+ 




[ah + ah) ~ ~2 Jq 


du 


(C.5) 


Next, we employ Holder’s inequality with the pair (p, g), where p, g > 1 and g = p/(p — 1), to 
force the term f^du outside the expectation. Then, we condition the second expectation 
on the cr-algebra and proceed as in (C.4) to arrive at 


E [L^] <S^E 


exp I poj Tudu 


/ 


E 


exp < q(jjai 4 ^ 


I 


* rz 

'Vu -f- du 
ot 


+ q 




f 


V,j du 


(C.6) 


All that we have left to do is to show that the supremum over t and 5t of each of the two 
expectations on the right-hand side of (C.6) is finite. However, one can easily deduce from 

(C.7) 


Lemma [3T] the following sufficient conditions: 


kd > k > qujpsv^ + ^ 


where 


A = q 


^oj{u - 1) + - l)(a?3 + ah) 


(C.8) 


Note that we used J2j=i ^ij = 1 in (C.8). On the other hand, the first assumption in (3.10) 
ensures that ^ 2 ( 0 ^) > 1- This, in turn, implies that q 2 {a)/{q 2 {a) — 1) > 1. Due to the second 
assumption in (3.10), we can find p > 1 so that 

aq2{a) 


— 2 > ^P> 


Ti, 


d 




(C.9) 


The quadratic equation in x below has roots of different signs, with positive root q 2 [c()- 
1 


-Jahh<^is + ah)^^ + 


otpsv^ + ^Ta^^^[ah + ah) “ 


X — A: = 0. 


However, the Holder pair satisfies p = q/{q — 1), so q < q 2 {o^). Hence, q lies in between the 
two roots of the quadratic, which implies that 


iT(a?3 + 0^4) q^ + apsv^ +^T[ah+ah) 
Rearranging terms in the above inequality, we obtain 


q — k <{). 


1 9 

k > qapsv^ + Q 


-a(a - 1) + -a‘^{q - l){ah + ah) 


(C.IO) 


(C.ll) 
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From (C.9) and (C.ll), employing a continuity argument similar to that used in the proof 


of Proposition 3.2, we deduce that the two conditions in (C.7) hold on an interval [a, 0 ^ 2 ), 
for some 0^2 > (a, which concludes the proof. The extension to the interval [ 1 , 0 ^ 2 ) follows 
from Jensen’s inequality, with r]^ — r]a^ Vcj G [1, a]. ■ 

In the event that ais and au are simultaneously zero, i.e., psv — Psd — 0, one can easily 
show that Proposition |3.3| still holds as long as 


1 


k > — l)T^ , 


> 


2ak 


^"Ak-a{a- l)r^2 


Appendix D. Proof of Proposition 


3.4 


We follow the argument of Proposition closely and condition on the a-algebra instead 
to deduce that 


rt rt 

Vudu-\-(jjpsv / 

/o Jo 

First of all, suppose that a = 1 and T > 0. If fc < Psv^^ then 


E[i?^]<S'oE exp I j Vudu + uips^J 


(D.l) 


1 , fujpsv^-k + u{uy 

-?—r- tT ' = 

0 ) 4 , 1 + viLO) \L0Psv^ — k — u{uj) ^ 

Hence, by a continuity argument, we can find ai> 1 such that for all uj G (1, cri), 

, ^ !—7-^ 1 , (UJPsvt — k^v[uj\ 

k < ujpsv^ - yw(w - 1)^ and T < —— log ( 

u{u) V 

If k = Psv^i then psv € (0,1] and 
2 


ujpsvi - k- u{uj) J ■ 


(D.2) 


lim 


TT / Ujpsv^ - k 

— — arctan 


1-2 ( Wlo — (lli — l)pL 

= lim —r arctan ( ^ 


0)4,1+ i>{oj) 


o;Tl+ k{L0) [2 
Furthermore, note that for all a; > 1, 

‘^Psv^ — 'J— 1) ^ < ojk — \/cj{oo — 1) k < k. 
Hence, we can find cii > 1 such that for all cj € (1, ai), 

ujpsv^ - - 1)^ < k < Lopsv^ + A/n;(w - 1)^ 

and 


■s/ui^ Ps 


= oo. 


T < 


[2 


TT f Upsvi - k 

— — arctan 


V 


(D.3) 

(D.4) 
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If A; > Psvii then we can find ai> 1 such that for all cj G (1, ai), 


k > IjJPsvi + a/w(cj - 1)^. 


(D.5) 


The conclusion follows from Proposition 3.2 in Cozma and Reisinger (2015b) and (D.2) - 
(D.5). Next, suppose that a> 1 and T < with T* defined in (3.16) - (3.19). 

If A: < apsvi — ^/a{a^^T) by a continuity argument, we can find ai > a such that for 
all (jj G (cr, CTi), 


k < ojpsv^ - and T < log ( 

iy{Lo) V 


Lopsv^ - k + u{Lij) 


^UPsv^ - k - viui) J ■ 

If k = apsvi — \/o({a — 1) then psv € (0,1] and for all lo > a, 

Lo - a < ^a{a - 1) copsv^ - \/oj{uJ - 1)^ < apsv^ - \/a{a - 1)^. 

Furthermore, note that 


(D.6) 


lim 


ujia+ P{oj) [2 


TT (UPsv^ - k 

— — arctan 


V H<^) 


= lim 




arctan 


Hence, we can find ai > a such that for all cj G (cr, cri). 


and 


T < 


TV ( bJpsvi - k 

— — arctan 



2 

V ^Psv^ ~ k ) 

(^Psv^ k 

UJ - 1)^ 

(D.7) 

- 

(D.8) 


[2 \ 

If apsv^ — — 1) C < k < apsvi + C? we can clearly find oii > a so that both 

(D.7) and ( |D.8 ) hold for all uo G (cr, ai). If A; = o^Psvi + then for all u > 

{a - uj)psv < - 1 ) - apsv^ + ^ 0(0 - 1 )^ < copsv^ + \/oj{u - 1 )^. 

Furthermore, note that 

f ^ psv^ k 


lim 


TV 

-arctan , ^, 

2 V 


r 

= hm , = oc. 

u;^a+ VyUj) 


Hence, we can find ai> a such that both ( |D.7 ) and (D.8) hold for all uo G (a, ai). Finally, 
if A: > apsv^ + (^, then we can find ai > a such that for all uj G (a, ai). 


k > Lopsvi + - 1)^. 


(D.9) 


The conclusion follows from Proposition 3.2 in Cozma and Reisinger (2015b) and (D.6) 
(D.9). The extension to the interval [l,(ai) follows from Jensen’s inequality. I 
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Appendix E. Proof of Proposition 3.5 


We follow the argument of Proposition |3 . 3| closely and condition on the cr-algebra instead 
to deduce that 


E [R^] <S^E 


exp 


pIv) - 


/ Vudu-\- ujpsv / Y v„, —7T^ du 
/o Jo 


St 


• (E.l) 


Suppose that T < with T* from (3.21) - (3.22). If < apsvi + \ \/a{a — 1) (^, then by 
a continuity argument, we can find a 2 > a such that for all cj G (cr, cr2). 


k < upsvi + - 1)C and T < --- ^ --. 

^ + Vw(cj - 1) ^ - A: 

If A: = OiPsv^ + \ since A:, ^ > 0 and for all ui > a, we have 


(E.2) 




Furthermore, note that 
4(A: - apsvi) _ 


and lim 


cr(cr-l )(^2 apsvi + cr(cr - 1)- A: ujpsvi + -l)^-k 


= r*. 


with T* from ( |3.22 ). Hence, we can find cr 2 > cr such that (E.2) holds for all uj G (cr,cr 2 ). 
Finally, if A; > apsv^ + ^ y/cr(cr — 1)(^, since 

lim 

uia+ 0J[UJ — 1 ) 4 "' 


with T* from (3.22), we can find 012 > a such that for all u € (a, 012 ), 


k > ujpsv^ + t yw(w - 1) ^ and T < 

2 6c;(6c; — 1)^^ 


(E.3) 


The conclusion follows from Lemma 3.1 and (E.2) - (E.3). The extension to the interval 


[1, cr 2 ) follows from Jensen’s inequality, with Pu = Vcj G [1, cr]. 


Appendix F. Proof of Proposition 


3.6 


The following auxiliary result proves the almost sure positivity of the foreign interest rate. 
Lemma F.l. Let k > (^ 0 )“^ and define the stopping time 

T,^ — inf > 0 : r/ < (FT) 
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If2kf9f > then 


lim P < r) =0. 


Proof. Define the function U : (0, oo) M by 


U{x)^x-r a^^{2kfef-ef). 


(F.2) 


(F.3) 


By Ito’s formula, we have 


E 


^U{rl)-E a{rl)~^^^'"\kfef-kfrl - psAf\lvsri)ds 

J 0 

+ -E (F.4) 

However, 

pT ^ ^ , o ^ 

ts<r^ ds < < oo, 


E 


[ ^ .2c2^1+2a. 

Jo 


so the stochastic integral on the right-hand side of (F.4) is a true martingale. Hence, 

cTAt^ 


E 

where 


{“(d) (F.5) 


“ = - «/)'• '> = -f/),c= ^(2kfDf - (}). (F.6) 


Employing Fubini’s theorem and Holder’s inequality in (F.5), we get 

cT 


E 


v(4aJ <i'(4)- 


c sup E 
nG[0,r] 


/o 


aE 




,1+c^ 


2(l + a) 


E 


( 4 ) 


s<r^ 


-(l+a) . 


-bE 




5 <Tk 


1 + a 


^S<T^ 


l+2a 
2(l + a) 


ds. 


(F.7) 


The moments of the square root process are uniformly bounded (Dereich et al. 2012) and 
the function / : [0, oo) M defined by 


1+20; 


f{x) — ax — bx^+^ — c sup E 2 ( 1 + 0 ;)^ 2 (i+a) 

nG[0,T] 


(F.8) 


is clearly bounded from below. Hence, we can find a constant C independent of k such that 


E 


U(rLr.) 


< c. 


(F.9) 


Since has continuous paths, we have = k ^ and U{rlA = Therefore, using (F.9) 


/ 'i - 
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and the fact that U is positive, we deduce that 


K'^¥{T,<T)=E\u{4jlr^<: 




u{4aJ ir.<Tl < E \u{4^j] < a (F.IO) 


Letting ^ oo in ( |F.10[ ) yields the conclusion. ■ 

The next two lemmas give moment bounds for the original and the discretized foreign 
interest rate processes. 

Lemma F.2. The process has uniformly bounded moments, i.e., 


E 


sup {r{y 

t&[0,T] 


< oo, Vp > 1. 


Proof. Fix any p > 1. From 

4 ^ 4 + [ 4du-ps4f[ \Jvuridu + if [ ffdwl. 

Jo Jo Jo 

Using the fact that 2^/\ab\ < \a\ + \b\ and Holder’s inequality, we deduce that 

{r{Y < 22(p-i) [4 + kfdfty + 2P-‘^\psf\Py ^ VudJj 

+ 2P-\2kf + \Pss\^fy{^j\idJj + 2^^P-^\P fJudWl 


(F.ll) 


(F.12) 


(F.13) 


Fix t G [0,r]. Using Holder’s inequality, we get 


sup {4y < 2 ^^P-^\rl + kfOfTy + 2P-^\p,f\PefTP-^ rvPdu 

se[o,t] Jo 

+ 2 P- 42 kf + M^fyTP-'^ [\4ydu+2^(p-'^kp sup 

Jo sG[0,t] 



From the Burkholder-Davis-Gundy inequality, we know that there exists a constant Cp > 0 
such that 


E 

sup 

r Jddwl 

p- 

<CpE 

If r f 

/ rLdu] 


. sG[0,t] 

Jo 



L \ ^0 / J 


1 


<^Cp + -Cj,TP-^E 


I Jo 


yiydu 


Taking expectations and employing Fubini’s theorem in (F.14), 
E 


sup {4y <2‘^^p-^y4 + kfefTy+2‘^p-4^fCp+2P-^\psf\pefTP sup e[<] 
se[o,t] J ue[o,T] 

+ ( 2 P- 42 kf + \p,f\^fyTP-^ + 2^P-^efCpTP-4 ['eI sup {4y du. 
^ ^ Jo Lsg[ 0 ,«] 
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Applying Gronwall’s inequality, we get 


E sup {r{y <(2^(P-^y4 + kf9fTY + 2^P-^efCp + 2P-^\psf\PefTP sup E[<]) 
LtG[0,r] J ^ nG[0,r] ^ 

X exp {2P-^2kf + IPsfl^fYTP + 2^P-^efCpTP] . (F.15) 

The conclusion follows from the boundedness of moments of v. ■ 

Lemma F.3. The process from (2.9) has uniformly hounded moments, i.e., 


sup E 
(5tG(o,7y) LtG[o,r] 


sup YlY 


< oo, Vp > 1, Mr] > 0. 


(F.16) 


Proof. Fix any p > 1 and p > 0. From (2.8), 


~f ^ f 


' t — • {) 


^ ^ I f I 

Since < |f/ |, following the argument of Lemma 
E 


ra + 


kfOft-kf [ rldu-psfCf [ ^Jvpldu + ^f [ ffdWY (F-17) 
Jo Jo Jo 


F.2 


we deduce that 


sup {fiY <2^^p-^Y4 + kfefTY+2^p-^efCp+2P-^\psf\pefTp sup Ey^^’] 
se[o,t] J ue[o,T] 

t 


+ (2P-^2kf + \p,f\^fYTP-^ + 2^P-^efCpTP-^) f E [{riY\ du. (F.18) 


Since sup„g[o,T] [^] ^ sup„g[o,r] V defined as in (^), and < sup^gjo,^] > 

applying Gronwall’s inequality, we get 


E 


sup (r/y < + + + sup Ey^^’]) 

tG[o,r] J ^ nG[o,r] ^ 

X exp {2P-^2kf + \Psf\^fYT^ + 2^P-^efCpTP] . (F.19) 


The conclusion follows from Proposition 3.4 in Cozma and Reisinger (2015a). 


Next, we use Lemma F.3 to prove the convergence of the L? difference between the two 
time continuous discretizations. 


Lemma F.4. The difference between r^ and r^ converges to zero with 6t, i.e., 


lim sup E 

tG[o,r] 


/-/ -/\2 

Yt - n) 


= 0. 


(F.20) 


Proof. Suppose that t € [t„,fn+i)- Since |r/ —r{| < \f( — f( \ and from (2.8), we can bound 
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the squared absolute difference from above as follows: 

{H - Hf < + 0.b\p,f\^fStVt^ + {kf + 0.b\p,f\Cf)6tfl + ^/y^|w/ - w/j) 

< Ak}0}{5tf + + {2kf + \psf\^f)\5tf{r(S + {w/ - Wl)\ 


Therefore, 


sup E\{r(-r{f] < Akjej{5tf + \psf\'^ff{5tf sup E 
te[o,r] L ^ te[o,r] 


+ {2kf + \psf\^f) {St)^ sup E {r{J +Aff5t sup E[r{J. (F.21) 

0<n<N *- 0<?T,<W 


Using Lemma F.3 as well as Proposition 3.4 in Cozma and Reisinger (2015a) concludes the 
proof. ■ 

The following lemma derives the strong mean square convergence of the stopped process. 

Lemma F.5. Let I > and define the stopping times 

Ti — inf {t > 0 : 'T’t > /} and r = T/^At/, (F.22) 


with defined in (F.l). Then the stopped proeess eonverges uniformly in , i.e., 

= 0 . 


lim E 
St^o 


sup 

te[o,T] 




(F.23) 


- I- , f f 

Proof. From (F.12) and (F.17), since |r^ ~ — rt ~ U I? have 


I / -f 
IUat Uat 


+ 


pt/\T rtAr / ! - /-\ 

-kf J {rl - fl)du -kf J {fl - ri)du + ^f j (y rl - ^ fljdWi 

if {\fii-fpjdWl-Psfif ' 


rtAr 

- psfif / 

JO 


fl - yrljdu- psfif 



. (F.24) 


Fix t G [0,T]. Squaring both sides and using Cauchy’s inequality, then taking expectations 
and using Doob’s martingale inequality and Fubini’s theorem, we get 


du 


E 

sup (Aat-Aat) 

< 7k)T [ E 


A 

du + 7k)T / E 

{fl-df 


_ sG [0,t] 

Jo 


Jo 



+ 


28«? f 


E 


) kn<r 


du + 


284? /’ 
^0 


E 


If/ _ f/1 

\' U ' u\ 


du 
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+ 7plfefT / E 

Jo 


Vu[yri 


'^u ) ln<r 


du + 7p‘ 


;ff,T f 


E 


ri\vu - Vu 


+ 7p: 


\,e,T f 


E 


Vu\rl-rl\lu<T 


du. 


du 
(F.25) 


On the other hand, we know that \rl — rl \ lu<r ^ ~ '^nArl 


(pfdx-\[HJj l«<r < (\4^- < n{rLr - rlArf. (F.26) 


Substituting back into (F.25) with (F.26), we arrive at the following inequality: 

rt 


E 


/ f -/ \2 
sup (Aat-Aat) 

_ sG [0,t] 


< {7k}T+ 28^^K + lplf^fTln) [ 

Jo 


E 


/ / \2 
sup [riAr-riAr) 

sG[0,n] 


du 


+ 7k}T^ sup E {rl-rlf + (28e^T + T^/) sup E 

ue[o,r] L ^ ue[o,T] 


If/ _ rf\ 

\' U ' u\ 


+ 7plf^‘jT‘^ sup E (rlf 
«e[o,T] 


sup E 
nG[0,r] 


I t7 | 2 ' 


(F.27) 


The convergence to zero of the last three terms on the right-hand side of (F.27) follows from 
Lemmas F.3 and F.4, and Proposition 3.5 in Cozma and Reisinger (2015a). The conclusion 
follows from a simple application of Gronwall’s inequality. ■ 


From Lemma [F.4| we know that in order to establish the strong mean square convergence 
of it suffices to prove this for r-^, since 

\r{-f{f <2\r{-f{\^+ 2\f{-f{f, Vt€[0,r]. (F.28) 

Lemma F.6. If2kf6f > then the process converges strongly in L?, i.e., 


lim sup E 

5t —>-0 t^[0,T] 


f V|2' 


r; - r: 


= 0 . 


(F.29) 


Proof. Fix K > (rg) I > uq, and recall the definition of the stopping time r from (F.22). 
Since r^ and are non-negative, 


sup E 

'CN 

1 

1_ 

< sup E 

1 / 1^ 11 
Vt-n\ lr<t 

-|- sup E 

1 / "/ 1 
Vt -n \ 

te[o,T] 


te[o,r] 


te\0,T] 



H<T 


< sup E 
te[o,T] 


{rlftrKT 


+ sup E 
te[o,T] 


(r/)^lr<r 


+ sup E 

tG[0,T] 


I / _ |2' 

I At At | 
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Since lr<T < + '^ ti<t and applying Cauchy’s inequality, we get 


sup E 

ll / 

In - 

-/121 
-N 1 

< < 

1 sup E 


1 

^ -b sup E 


t€lO,T] 




1 t€lO,T] 


tG[0,T] 






sup E 


^ + sup E 






1 te[o,T] 


te[o,r] 



■ (t- < T) 

{n<T)‘- 


+ sup E 
te[o,T] 


/ _ -/ |2' 


(F.30) 


On the other hand, using Markov’s inequality, we obtain an upper bound 

1 


^ {tI < T) < F ( sup Vf > I \ < y E sup Vt . (F-31) 

\t€lO,T] ) I LtG[0,T] _ 

However, the expectation on the right-hand side is clearly finite by the Burkholder-Davis- 


Gundy inequality. Taking the limit as St ^ 0 in (F.30) and employing Lemmas F.l to F.3 
and |F.5[ since k and I can be made arbitrarily large, leads to the conclusion. ■ 


Appendix G. Proof of Proposition 


3.7 


Note that Vt € and Vj € {2,3,4}, since V is piecewise constant, 


n fir ^ fir - Wi f^Wl - W/ 


z=0 


VudWi^sJVt 


i — i. 






^Jf+l _ t tu/i 


6t 


(wi - o 


For convenience, Vt G [tn^tn^i) and Vj G {2,3,4}, we define 

t tr). ^ 


^’= St 


-in, - wn - - wi). 


(G.l) 


(G.2) 


Substituting back into (3.25) with (G.l) and (G.2), we obtain 


rt -1 rt ! - ^ /- 

Xt = xo+ {ri-ri--Vu)du+ Vu dW^ + Y, aij Vt Z{. 

^ 0 ^ 0 j 2 


(G.3) 


The absolute difference between the original and the discretized log-processes is thus 




[ {ri-rt)du- [ {rl-rl)du-j^ [ {vu-Vu)du 
Jo Jo ^ Jo 


rt , ^ ri 

+ iWt-yVu)dw^-Y «ii yv zi 

Jo ^ ^ i=2 


(G.4) 
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Squaring both sides of (G.4), applying the Cauchy-Schwarz inequality, taking the supremum 
over all t G [0,T], and then using Cauchy’s inequality for all Riemann integrals leads to 

rT pT 7 rT 


sup < 7T 

te[o,T] 


^ {rt-fiYdu + 7T (rl-rifdu+^-T {v^-V^fdu 


+ 7 sup 

tG[0,T] 


^-\/Vu dW^ 


+ 7y^a^^- sup Vt sup |^/|^. (G.5) 

teio^n telo,T] 


We used the fact that sup^^jQj^j Vt < sup^^jQj-] Vt^ with V defined as in (2.6). Taking expec¬ 
tations and employing Fubini’s theorem, Holder’s inequality, Doob’s martingale inequality 
and the Ito isometry, we derive 


E 


1 |2 
sup \xt-Xt\ 

< 7T'^ sup E 

]rf 

-nl 

+ 7T‘^ sup E 

[l / -/|2l 

rt -ri\ 

Itelo.T] J 

te[o,T] 



tG[0,T] 



7 9 

+ - sup E 
4 te[o,T] 


[l t7|21 

\vt - Vt\ 

-|- 28r sup E 

1 

_1 


tG[0,T] 



+ 7 afj E 
i=2 


sup Vf 
Lte[o,T] 


1/2 


E 


sup I Zj 
te[o,r] 


1/2 


(G.6) 


The convergence as 5/ —> 0 of the first four terms on the right-hand side of (G.6) follows from 
Propositionand Proposition 3.5 in Gozma and Reisinger (2015a). Next, integrating the 
time continuous auxiliary variance process defined in (2.5) leads to 

rt ft I — 

(G.7) 


vt = vo + k f {e-Vu)du + ^[ JVudWt 
Jo Jo 


However, V = max{0, v} < |i)| and, using Cauchy’s inequality, Fubini’s theorem and Doob’s 
inequality, we find an upper bound 


E 


sup V{ 

LtG[0,T] 


^ < 3(^o + A;0T)V3fc2r2 sup ^[V^]+l2i^T sup E[Vt]. (G.8) 

te[o,T] te[o,T] 


The uniform boundedness of the second moment of the FTE discretization for the variance 


as ^ 0 follows from Proposition 3.4 in Gozma and Reisinger (2015a). Finally, employing 
the definition in (G.2), we bound the term inside the last expectation in (G.6) from above. 

t-tr 


sup sup sup 

te[0,T] 0<n<N tn<t<tn-\-l 

< 8 sup sup 

0<n<N tn<t<tn-\-l 


6t 


- w'i) - W - O 

4 


t — t 
St 


+\^i -n. 


i |4 
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<16 sup sup \Wl — 

0<n<N tn<t<tn-\-l 

<16 ^y \wl-wl i \\ 

te[o,r] '■' 


1 |4 


(G.9) 


However, moments of the Euler modulus of continuity of a Brownian motion converge to 0 
as ^ 0 ( Fischer and Nap^^|2009 ), which concludes the proof. ■ 


Appendix H. Proof of Proposition 3.8 


For fixed, positive numbers e and 7 such that log (1 + e) > 7 , define the set 

Se ,7 = {x G M I 3^ G M : |x — ^| < 7 and |e^ — e^| > e}. (H.l) 

However, since the exponential function is strictly increasing, 

xGHe,<y 47 3^G (x-7,x + 7): - l) > 6 47 e^(e^-l) > e. 

Hence, 

= (a( 6 , 7 ),+ 00 ), where a(e, 7 ) = log (^ 7 ^) > 0. (H.2) 

We have the following string of inclusions of events. 


sup I > e > C < sup —Xt|>7>U< sup \xt — Xt\<^^ sup |e^^ —e^^|>< 

telo.T] J Ug[o,t] J UG[o,r] te[o,T] 


c sup Xt-Xt >7 > U < G [0,r] : xt G 
.tG[o,r] 


C< sup \xt — Xt\>j>U< sup xt>a(e,7 ) 

UG[o,r] J UG[o,r] 

In terms of probabilities of events, the previous inclusion becomes: 


(H.3) 


)+p( sup xt>a{e,'y)]. (H.4) 

/ \t&[0,T] 


sup pt — S'*! > e I < P I sup \xt — Xt\ >7 
te[o,T] J \te[o,T] 

The convergence in probability of the log-process is a consequence of Proposition |3.7| and 
Markov’s inequality. Therefore, all that we have left to prove is that the second probability 


on the right-hand side of (H.4) can be made arbitrarily small. However, if we fix e > 0 and 


vary 7 , then lim.y_j.o a(e, 7 ) = 00 . A simple application of Markov’s inequality leads to 


sup Xt > a(e, 7 ) ) < 
tG[0,T] 


sup Ixtl > a(e, 7 ) ) < 3 . 

te[o,T] / a(e,7) 


E 


sup \Xt\ 
te[o,T] 


(H.5) 
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On the other hand, nsing Jensen’s ineqnality and Doob’s martingale ineqnality, 


E 


snp I a;* I 

LtG[0,T] 


< |a:o| + sup E \rf\ + T sup E [r/] + 


te[o,T] 


te[o,T] 


sup 

te[o,T] 


E [vt] + 2VT sup E [vt] 


te[o,T] 


However, the right-hand side is finite because the moments of the square root process are 
bounded and from Proposition F.2, which concludes the proof. ■ 


Appendix I. Proof of Theorem 


3.9 


Fix e > 0 and define the event A = > e|. Since S and S are non-negative, 


sup E 

1 

P 

< sup E 


+ sup E 


te[o,r] 


te[o,T] 


te[o,r] 



< -h sup E[S^1a]+ sup E[S^1a]- (FI) 

tG[o,r] telo,T] 

Let a < (jj < min{(ai,(a 2 } and apply Holder’s inequality to the two expectations on the 
right-hand side of p.l| ) with the pair (p^q) = (^, Hence, 


sup E 

telo,T] 


pt - St 


< e“ + 


sup E[S'^]‘^+ sup E[S'^]‘^1 sup P — ^tl > e 
te[Q,T] te[o,T] J te[o,T] ^ 


The convergence of S in probability is a consequence of Proposition 3.8 Finally, employing 


Propositions |3.2| and |3.3| and then taking e sufficiently small concludes the proof. 


52 















